Despite being large molecules, DNA and proteins are made with pieces of only a few types

DNA is made of four bases. We can represent it with four letters

AGCTTTTCATTCTGACTGCAACGGGCAATATGTCTCTGTGTGGATTAAAAAAAGAGTGTCTGATAGCAGC

Proteins have 20 aminoacids. We can represent them with 20 letters.

MRVLKFGGTSVANAERFLRVADILESNARQGQVATVLSAPAKITNHLVAMIEKTISGQDALPNISDAERI

In molecular biology we often work with sequences

• DNA sequences use 4 letters to represent the nucleotides in one of the two strands
• Protein sequences use 20 letters to represent the amino-acids, from amino to carboxyl terminal
• Other sequences are sometime used:
  • RNA,
  • DNA with ambiguous nucleotides,
  • amino-acid sequences with stop codons

The main reason why computing is useful for molecular biology

• DNA is discrete data
• Either “A”, “C”, “G” or “T”
• No middle-values
• All other things we measure are in a range
  • Like temperature, concentration, expression
• We can forget all the chemistry and work with symbols

There are several ways to store DNA or protein data

Most of the times they are stored in FASTA format

FASTA files are text files, with some rules

Microsoft Word files (doc or docx) are NOT text files

You should never use Microsoft Word to store sequences

We use this genome in classes because it is a small example

Candidatus Carsonella ruddii is a obligate symbiont of Pachpsylla venusta.

There is much controversy over whether this is a living cell or simply an organelle as it is missing genes needed for living independently.

Published as “The 160-kilobase genome of the bacterial endosymbiont Carsonella.” https://www.ncbi.nlm.nih.gov/pubmed/17038615

NCBI stores all public biological sequences at
https://www.ncbi.nlm.nih.gov/nuccore

Anybody can upload sequences, and they may be wrong

NCBI has a curation process to validate the sequences

If a sequence is good enough to be a reference, then it is stored in the RefSeq collection https://www.ncbi.nlm.nih.gov/refseq

Accession numbers are the best way to identify a biological sequence

• Candidatus Carsonella ruddii PV DNA has accession AP009180
• Escherichia coli str. K-12 substr. MG1655 has accession NC_000913

Different sequences can have the same name, but never the same accession

NCBI search box looks for patterns everywhere, not only on the title

and many sequences can have the same name

To be sure of using the correct sequence, go to NCBI Taxonomy

https://www.ncbi.nlm.nih.gov/taxonomy

Try with “Escherichia coli”. There are many

To handle sequence data in R, we use the seqinr library

You have to install it once.

install.packages("seqinr")

Then you have to load it on every session

library(seqinr)

read.fasta(file, seqtype = "DNA", ...)

Inputs

file

The name of the file which the sequences in FASTA format are to be read from

seqtype

the nature of the sequence: DNA or AA, defaulting to DNA

Output

A list of vectors of chars. Each element is a sequence object.

• Some organisms have been sequenced
• We can read the sequences from FASTA files
• read.fasta returns a list of vectors of characters
• We can create functions to do things with the sequence
• For example we do statistics of the sequence

Statistics is a way to tell a story that makes sense of the data

In genomics, we look for biological sense

That story can be global: about the complete genome

Or can be local: about some region of the genome

We will start with global properties

Measuring the melting temperature of the DNA double helix using spectrophotometry

• The absorbance of DNA at a wavelength of 260nm increases when DNA separates into two single strands at melting temperature

If the DNA has been sequenced then the GC-content can be accurately calculated by simple arithmetic.

GC-content percentage is calculated as \[\frac{G+C}{A+T+G+C}\]

Write the code for this function:

genome_gc_content <- function(genome) {
    # count letters
    return(answer)
}

Then use it to find the GC content of E.coli

Is the GC content uniform through all genome?

Do all genes have the same GC content?

What function do we need to answer this question?

What is the name, the inputs, and the output?

Write the code for this function:

genes_gc_content <- function(genes) {
    # calculate GC content of each gene
    return(answer)
}

What is the data structure of answer?

Calculate the GC content for only part of the genome?

Instead of all the genome, we only look through a window

The result should depend on:

• the genomic sequence
• the position of the window
• the size of the window

The new function looks like this

gc_content <- function(sequence, position, size) {
    # code here
    return(answer)
}

What should be the code?

We want to evaluate gc_content on different positions of the genome

positions <- seq(from=1, to=length(genome), by= window_size)

Now we apply gc_content to each element of positions?

for(i in 1:length(positions)) {
    gc[i] <- gc_content(genome, positions[i], window_size)
}