!curl -Lo conda_installer.py https://raw.githubusercontent.com/deepchem/deepchem/master/scripts/colab_install.py
import conda_installer
conda_installer.install()
!/root/miniconda/bin/conda info -e

!pip install --pre deepchem
import deepchem
deepchem.__version__

"""
Explanation: Tutorial Part 21: Introduction to Bioinformatics
So far in this tutorial, we've primarily worked on the problems of cheminformatics. We've been interested in seeing how we can use the techniques of machine learning to make predictions about the properties of molecules. In this tutorial, we're going to shift a bit and see how we can use classical computer science techniques and machine learning to tackle problems in bioinformatics.
For this, we're going to use the venerable biopython library to do some basic bioinformatics. A lot of the material in this notebook is adapted from the extensive official [Biopython tutorial]http://biopython.org/DIST/docs/tutorial/Tutorial.html). We strongly recommend checking out the official tutorial after you work through this notebook!
Colab
This tutorial and the rest in this sequence are designed to be done in Google colab. If you'd like to open this notebook in colab, you can use the following link.

Setup
To run DeepChem within Colab, you'll need to run the following cell of installation commands. This will take about 5 minutes to run to completion and install your environment.
End of explanation
"""


!pip install biopython

import Bio
Bio.__version__

from Bio.Seq import Seq
my_seq = Seq("AGTACACATTG")
my_seq

my_seq.complement()

my_seq.reverse_complement()

"""
Explanation: We'll use pip to install biopython
End of explanation
"""


!wget https://raw.githubusercontent.com/biopython/biopython/master/Doc/examples/ls_orchid.fasta

"""
Explanation: Parsing Sequence Records
We're going to download a sample fasta file from the Biopython tutorial to use in some exercises. This file is a set of hits for a sequence (of lady slipper orcid genes).
End of explanation
"""


from Bio import SeqIO

for seq_record in SeqIO.parse('ls_orchid.fasta', 'fasta'):
    print(seq_record.id)
    print(repr(seq_record.seq))
    print(len(seq_record))

"""
Explanation: Let's take a look at what the contents of this file look like:
End of explanation
"""


from Bio.Seq import Seq
from Bio.Alphabet import IUPAC
my_seq = Seq("ACAGTAGAC", IUPAC.unambiguous_dna)
my_seq

my_seq.alphabet

"""
Explanation: Sequence Objects
A large part of the biopython infrastructure deals with tools for handlings sequences. These could be DNA sequences, RNA sequences, amino acid sequences or even more exotic constructs. To tell biopython what type of sequence it's dealing with, you can specify the alphabet explicitly.
End of explanation
"""


my_prot = Seq("AAAAA", IUPAC.protein) # Alanine pentapeptide
my_prot

my_prot.alphabet

"""
Explanation: If we want to code a protein sequence, we can do that just as easily.
End of explanation
"""


print(len(my_prot))

my_prot[0]

"""
Explanation: We can take the length of sequences and index into them like strings.
End of explanation
"""


my_prot[0:3]

"""
Explanation: You can also use slice notation on sequences to get subsequences.
End of explanation
"""


my_prot + my_prot

"""
Explanation: You can concatenate sequences if they have the same type so this works.
End of explanation
"""


my_prot + my_seq

"""
Explanation: But this fails
End of explanation
"""


from Bio.Seq import Seq
from Bio.Alphabet import IUPAC

coding_dna = Seq("ATGATCTCGTAA", IUPAC.unambiguous_dna)
coding_dna

template_dna = coding_dna.reverse_complement()
template_dna

"""
Explanation: Transcription
Transcription is the process by which a DNA sequence is converted into messenger RNA. Remember that this is part of the "central dogma" of biology in which DNA engenders messenger RNA which engenders proteins. Here's a nice representation of this cycle borrowed from a Khan academy lesson.
<img src="https://cdn.kastatic.org/ka-perseus-images/20ce29384b2e7ff0cdea72acaa5b1dbd7287ab00.png">
Note from the image above that DNA has two strands. The top strand is typically called the coding strand, and the bottom the template strand. The template strand is used for the actual transcription process of conversion into messenger RNA, but in bioinformatics, it's more common to work with the coding strand. Let's now see how we can execute a transcription computationally using Biopython.
End of explanation
"""


messenger_rna = coding_dna.transcribe()
messenger_rna

"""
Explanation: Note that these sequences match those in the image below. You might be confused about why the template_dna sequence is shown reversed. The reason is that by convention, the template strand is read in the reverse direction.
Let's now see how we can transcribe our coding_dna strand into messenger RNA. This will only swap 'T' for 'U' and change the alphabet for our object.
End of explanation
"""


messenger_rna.back_transcribe()

"""
Explanation: We can also perform a "back-transcription" to recover the original coding strand from the messenger RNA.
End of explanation
"""


coding_dna.translate()

"""
Explanation: Translation
Translation is the next step in the process, whereby a messenger RNA is transformed into a protein sequence. Here's a beautiful diagram from Wikipedia that lays out the basics of this process.
<img src="https://upload.wikimedia.org/wikipedia/commons/thumb/b/b1/Ribosome_mRNA_translation_en.svg/1000px-Ribosome_mRNA_translation_en.svg.png">
Note how 3 nucleotides at a time correspond to one new amino acid added to the growing protein chain. A set of 3 nucleotides which codes for a given amino acid is called a "codon." We can use the translate() method on the messenger rna to perform this transformation in code.
messenger_rna.translate()
The translation can also be performed directly from the coding sequence DNA
End of explanation
"""


coding_dna = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG", IUPAC.unambiguous_dna)
coding_dna.translate()

"""
Explanation: Let's now consider a longer genetic sequence that has some more interesting structure for us to look at.
End of explanation
"""


coding_dna.translate(to_stop=True)

"""
Explanation: In both of the sequences above, '*' represents the stop codon. A stop codon is a sequence of 3 nucleotides that turns off the protein machinery. In DNA, the stop codons are 'TGA', 'TAA', 'TAG'. Note that this latest sequence has multiple stop codons. It's possible to run the machinery up to the first stop codon and pause too.
End of explanation
"""


from Bio.Alphabet import generic_dna

gene = Seq("GTGAAAAAGATGCAATCTATCGTACTCGCACTTTCCCTGGTTCTGGTCGCTCCCATGGCA" + \
           "GCACAGGCTGCGGAAATTACGTTAGTCCCGTCAGTAAAATTACAGATAGGCGATCGTGAT" + \
           "AATCGTGGCTATTACTGGGATGGAGGTCACTGGCGCGACCACGGCTGGTGGAAACAACAT" + \
           "TATGAATGGCGAGGCAATCGCTGGCACCTACACGGACCGCCGCCACCGCCGCGCCACCAT" + \
           "AAGAAAGCTCCTCATGATCATCACGGCGGTCATGGTCCAGGCAAACATCACCGCTAA",
           generic_dna)
# We specify a "table" to use a different translation table for bacterial proteins
gene.translate(table="Bacterial")

gene.translate(table="Bacterial", to_stop=True)

"""
Explanation: We're going to introduce a bit of terminology here. A complete coding sequence CDS is a nucleotide sequence of messenger RNA which is made of a whole number of codons (that is, the length of the sequence is a multiple of 3), starts with a "start codon" and ends with a "stop codon". A start codon is basically the opposite of a stop codon and is mostly commonly the sequence "AUG", but can be different (especially if you're dealing with something like bacterial DNA).
Let's see how we can translate a complete CDS of bacterial messenger RNA.
End of explanation
"""


from Bio.SeqRecord import SeqRecord
help(SeqRecord)

"""
Explanation: Handling Annotated Sequences
Sometimes it will be useful for us to be able to handle annotated sequences where there's richer annotations, as in GenBank or EMBL files. For these purposes, we'll want to use the SeqRecord class.
End of explanation
"""


from Bio.SeqRecord import SeqRecord

simple_seq = Seq("GATC")
simple_seq_r = SeqRecord(simple_seq)

simple_seq_r.id = "AC12345"
simple_seq_r.description = "Made up sequence"
print(simple_seq_r.id)
print(simple_seq_r.description)

"""
Explanation: Let's write a bit of code involving SeqRecord and see how it comes out looking.
End of explanation
"""


!wget https://raw.githubusercontent.com/biopython/biopython/master/Tests/GenBank/NC_005816.fna

from Bio import SeqIO

record = SeqIO.read("NC_005816.fna", "fasta")
record

"""
Explanation: Let's now see how we can use SeqRecord to parse a large fasta file. We'll pull down a file hosted on the biopython site.
End of explanation
"""


record.id

record.name

record.description

"""
Explanation: Note how there's a number of annotations attached to the SeqRecord object!
Let's take a closer look.
End of explanation
"""


!wget https://raw.githubusercontent.com/biopython/biopython/master/Tests/GenBank/NC_005816.gb

from Bio import SeqIO

record = SeqIO.read("NC_005816.gb", "genbank")
record

"""
Explanation: Let's now look at the same sequence, but downloaded from GenBank. We'll download the hosted file from the biopython tutorial website as before.
End of explanation
"""


