Genome assembly preliminary results
1 / 62

Genome Assembly Preliminary Results - PowerPoint PPT Presentation

  • Uploaded on

Genome Assembly Preliminary Results. Jeri Dilts Suzanna Kim Hema Nagrajan Deepak Purushotham Ambily Sivadas Amit Rupani Leo Wu 02/01/2012. Outline. Data Pre-Processing Formats and Conversion PRINSEQ Data statistics Error Correction CoRAL Assembler results Newbler 2.6 MIRA3

I am the owner, or an agent authorized to act on behalf of the owner, of the copyrighted work described.
Download Presentation

PowerPoint Slideshow about 'Genome Assembly Preliminary Results' - ownah

An Image/Link below is provided (as is) to download presentation

Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author.While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server.

- - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - -
Presentation Transcript
Genome assembly preliminary results

Genome AssemblyPreliminary Results

Jeri Dilts

Suzanna Kim


Deepak Purushotham



Leo Wu



  • Data Pre-Processing

    • Formats and Conversion


    • Data statistics

  • Error Correction

    • CoRAL

  • Assembler results

    • Newbler 2.6

    • MIRA3

    • ABySS

    • Velvet

    • PCAP-454

  • Results

  • Lab Exercise

  • Outline1

    • Data Pre-Processing

      • Formats and Conversion

      • PRINSEQ

      • Data statistics

  • Error Correction

    • CoRAL

  • Assembler results

    • Newbler 2.6

    • MIRA3

    • ABySS

    • Velvet

    • PCAP-454

  • Results

  • Lab Exercise

  • What are sff files
    What are sff files? 

    • Sff files are Roche's "Standard Flowgram Format" files, containing the sequence data produced from a 454 run.

    • The sff files contains a Manifest header at the start describing the contents and flow intensity signal values for each base in each read.

    • They are in binary format, so need to be converted to text format, such as a fastq/fasta file using sff2fastq , ssf_extract , sffinfo programs. The Sequence Read Archive request that these .sff files be uploaded, to obtain accession number for publications.  

    Fastq fasta quality
    Fastq = Fasta + Quality

    • No standard file extension: but .fq .fastq and .txt are commonly used

      • 4 lines per sequence

      • Line 1 begins with the @ character, a sequence ID, and an optional description


      • Line 2 is the sequence letters


      • Line 3 begins with the + character, followed by the same sequence ID, and another optional description

       +Optional Description

      • Line 4 encodes quality values for the sequence letters in line 2


    File tools
    File Tools

    • There are a large number of sff tool converters available now. We list a few here:

    • sffinfo --> extract fasta, quality and flowgrams as text from .sff files.

      • sffinfo -seq Axolotl_reads.sff > Axolotl.fna

      • sffinfo -qual Axolotl_reads.sff > Axolotl.qual

    • sff2fastq --> extracts read information from sff and outputs the sequence and quality scores in fastq

      • sff2fastq -o  Axolotl.fastq  Axolotl_reads.sff

    • sfffile --> join sff files; extract part of sff file by MIDs, read names or random reads; or trim reads in user-defined ways.

    Quality control
    Quality Control

    WHY ?

    • Saves time, effort and money


    • Number and Length of sequences

    • Base qualities

    • Ambiguous bases

    • Sequence duplications 

    Data pre processing prinseq
    Data pre-processing - PRINSEQ

    Generating trimmed reads in fastq
    Generating trimmed reads in fastq

    • $ -fastq M19107.fastq -out_good stdout  -log M19107.txt -trim_qual_left 20 -trim_qual_right 20 -trim_tail_left 5 -trim_tail_right 5 -trim_ns_left 2 -trim_ns_right 2 -min_len 65 -min_qual_mean 20 -ns_max_n 4 | gzip -9 > M19107.fastq.gz

    Read a fastq file containing quality data and write data passing all Filters to standard out (terminal). The trimmed sequences are gzipped to a new file.

    -fastq:Fastq file containing sequence and quality data. 

    -out_goodstdout: This will write all data passing the filter to the stdout (terminal)

    -log:logfile to keep track of parameters, error etc.

    -min_len: Filter out sequences lower than this length.

    -ns_max_n: Filter sequences with more than the specified Ns. We tried with 2/3 and 1.

    -min_gc: Filter sequences with GC content less than min_gc. 

    -max_gc: Filter sequences with GC content greater than max_gc

     -min_qual_mean: Filter sequences with mean quality scores below the specified level. Most published thresholds varied between 15 and 25. We used 20. 

    Generating graphs report from trimmed reads
    Generating graphs report from trimmed reads -fastq M19107_filtered.fastq -graph_data -verbose -out_good null -out_bad null

    Reads a filtered fastq file and graph data file  to generate graphs showing the distribution of length, base quality, GC Content, Occurance of N, Poly A/T Tails, Sequence Duplication, Sequence Complexity and Dinucleotide odds ratios.

    -fastq:Fastq file containing sequence and quality data. 

    -graph_data: File containing graph data to generate graphs report 

    -verbose: prints status and info messages during processing

    -out_bad null & -out_good null: This will NOT create any output file other than the graphics file.


    • Data Pre-Processing

      • Formats and Conversion

      • PRINSEQ

      • Data statistics

  • Error Correction

    • CoRAL

  • Assembler results

    • Newbler 2.6

    • MIRA3

    • ABySS

    • Velvet

    • PCAP-454

  • Results

  • Lab Exercise

  • Error correction
    Error Correction


    • Sequencing errors pose the biggest challenge

    • Computational efficiency of assemblers improves 

    • Lot of redundant data - take advantage of it

    • Ensures high data usage in assembly

    Coral v1 3
    Coral v1.3

    • CORrection with ALignment

    • Corrects sequencing errors in short-read high throughput data

    • Key strategy - Multiple Alignment using redundant read data

    • Similar reads are all updated according to the alignment based on scoring of quality reads --score is based on number of reads that align at a position/ number of total reads at position

    Coral command line
    Coral Command Line

    • coral -f[q, s] Input File -o Output File 

    • accepts input files FASTA, FASTQ, and Solexa FASTA

    • -k (length of k-mer) >= log4 (length of genome), default 21

    • -e (maximum expected error rate), default 0.07

    • -454 (sets optimal settings in gap penalty, mismatch penalty, and reward for matching


    • Data Pre-Processing

      • Formats and Conversion

      • PRINSEQ

      • Data statistics

  • Error Correction

    • CoRAL

  • Assembler results

    • Newbler 2.6

    • MIRA3

    • ABySS

    • Velvet

    • PCAP-454

  • Results

  • Lab Exercise

  • Bler algorithm

    bler Algorithm

    How does it work??

    What is newbler
    What is newbler??

    • Roche's “GS De Novo Assembler” (where “GS” = “Genome Sequencer”)

    • Designed to assemble reads from the Roche 454 sequencer.

    • Accepts:

    • 454 Flx Standard reads, and

    • 454 Titanium reads.

    • Single and paired-end reads.

    • Optionally can include Sanger reads.

    • Runs on Linux, and has 32 bit and 64 bit versions.

    • Has Command-line and Java-based GUI interface.

    Olc algorithm a quick recap
    OLC algorithm -A quick recap

    • Overlap Layout Consensus (OLC) is a method used for de novo genome assemblies.

    • OLC requires three steps:

    • 1) overlap, 2) layout, and 3) consensus.

    • The overlap stage computes and builds the basic assembly graph.

    • The layout stage compresses the graph, and the consensus stage determines the genome sequence based on the graph generated in the previous two steps.


    In the overlap step, the sequence of each read is compared to that of every other read, in both the forward and reverse complement orientations.

    As such, the overlap computation step is a very time intensive step – especially if the set of reads is very large.

    Overlap criteria
    Overlap criteria

    • The OLC overlap criteria result in two types of overlaps: true overlaps (Figure 1a) and repeat overlaps (Figure 1b).

    • For example, in Figure 1b, an overlap occurs between reads S and T, due to the orange repeat section, not because reads S and T truly overlap one locus in the genome, as in Figure 1a.

    How does it go
    How does it go?

    • Although we would like exclude repeat overlaps, we must construct the assembly graph using both types of overlaps, as true and repeat overlaps cannot be distinguished individually.

    • In the assembly graph, the nodes represent actual reads, and edges represent OLC-quality overlaps between these reads (Figure 2).

    • Thus, the genome assembly becomes equivalent to finding a path through the graph that visits each node exactly once (i.e., a Hamiltonian path).


    • Finding a path through the OLC graph is not a trivial task. Imagine you have a graph of millions of nodes and edges.

    • Identifying a path that visits every node exactly once would be extremely difficult, even for a powerful computer.

    • In order to find such a path, you would have to start at some node and proceed to other connected nodes.

    • If you find that you visited a node more than once, you must backtrack, adjust the path, and test this new path.

    • So the larger the graph is, the more options you must test.


    • In order to decrease the size of the graph, the OLC assembly graph is simplified in the layout stage, where segments of the assembly graph are compressed into contigs, which are collections of reads that clearly overlap each other and refer to the same overall sequence.

    • Thus in the overlap graph, a contig would be a subgraph, or a group of nodes, with many connections among each other, as they all overlap with each other and refer to the same sequence

    Unique and repeat contigs
    Unique and repeat contigs

    • There are two classes of contigs, unique contigs and repeat contigs.

    • Unique contigs are composed of reads that can be unambiguously assembled.

    • Generally, these reads only overlap in one way and do not contain repeats within the genome. Essentially, unique contigs are contigs not flagged as repeat contigs.

    Contig assignment
    Contig Assignment

    Contig assignment1
    Contig assignment

    • Repeat contigs are contigs with an abnormally high read coverage or connected to an abnormally large number of other contigs.

    • Additionally, this repeat contig is different from other contigs because it has such high coverage.

    • The high coverage of repeat contigs allows algorithms to identify them through statistics that compare the coverage of each contig.

    • Contigswith too much coverage are most likely due to over-collapsing of repeats and are flagged as repeat contigs, to be used later in the layout stage


    • In the final stage of the OLC method, we derive the consensus sequence.

    Genome assembly preliminary results

    Identifying Unique DNA Stretches

    Repetitive DNA unitig

    Unique DNA unitig

    Arrival Intervals

    Discriminator Statistic is log-odds ratio of probability unitig is unique DNA versus 2-copy DNA.




    Dist. For Unique

    Dist. For Repetitive

    Definitely Repetitive

    Don’t Know

    Definitely Unique

    How to run newbler
    How to run Newbler?


    • The simplest command to run Newbler is:

    • runAssembly [options] reads.sff

    • • Which creates an the assembly in an output directory called:

    • P_yyyy_mm_dd_hh_min_sec_runAssembly

    • where P_ = Project, followed by date and time

    • There are a large number of optional parameters available for controlling and refining the assembly

    Lets look at a newbler run
    Lets look at a Newbler run..

    Overlap between two sequences
    Overlap between two sequences

    overlap (19 bases)

    overhang (6 bases)




    % identity = 18/19 % = 94.7%

    • overlap - region of similarity between regions

    • overhang - un-aligned ends of the sequences

    • The assembler screens merges based on:

    • length of overlap

    • % identity in overlap region

    • maximum overhang size.

    Newbler command line
    Newbler Command Line

    • runAssembly [options] inputFile

    • runMapping[options] referenceGenomeinputFile

    • accepts input files -- FASTA, FASTQ, .SFF

    • -o Name output directory

    • -qo Generate quick output for assembly (can decrease accuracy)  

    M19107 newbler results
    M19107 - Newbler Results

    M19501 newbler results
    M19501- Newbler Results

    Mira3 multifunctional inertial reference assembly v3 4 0
    MIRA3Multifunctional Inertial Reference Assembly v3.4.0

    • Started in 1997 as a PhD project at the German Cancer Research Center by BastienChevreux and in 2007 became open source

    • MIRA 3 is able to perform true hybrid de-novo assemblies using reads gathered through Sanger, 454, Solexa, IonTorrent or PacBio sequencing technologies.

      • can also perform regular (non-hybrid) de-novo assemblies using 454 data

    • Overlap layout consensus algorithm

    Mira3 command line arguments 1 2 extracting sff
    MIRA3:Command Line Arguments (1/2)"Extracting SFF"

    Convert the SFFs named M19107.sff, M19107.sff and M19107.sff

    The parameters of sff_extract

    -Q extract to FASTQ

    -s give the FASTQ file a name we chose 

    -x give the XML file with vector clipping information a name we chose

    Mira3 command line arguments 2 2 begin assembly
    MIRA3:Command Line Arguments (2/2)"Begin Assembly"

    Parameters--project (for naming your assembly project)--job (perhaps to change the quality of the assembly to 'draft')

    >& creates/outputs a file named log_assembly.txt to observe assembly progress

    Mira3 data manipulation integrity
    MIRA3 Data Manipulation Integrity

    Our objective is to produce the most accurate representation of a genome.

    When software tools produce better results, it doesn't necessarily indicate that the genome's representation is more accurate. (and vice versa)

    This makes it tricky to determine the proper tools to use.

    Could be detrimental to scientific integrity.

    Take precautions.

    Mira3 data analysis
    MIRA3 Data Analysis

    Data Quality (ideally) Increases

    Mira3 further work
    MIRA3 Further Work.....

    • Need to look at more assembly parameters in reference manual

    • Finish scripts that calculate Min. Contig Length and Avg. Contig Length

    • M19501 32bit error

      Finish running MIRA3 on all genomes


    • ABySS stands for Assembly ByShort Sequences.

    • It is a De Novo  sequence assembler designed for short reads and large genomes.

    • Single Processor version: Genomes up to 100Mbp in size  Parallel version: capable of assembling mammalian sized genomes

    • Capable of performing assemblies for both single end reads and paired end reads

    • The output of the ABySS is set of contigs assembled from the short reads

    Abyss commands for assembly
    ABySS commands for assembly

    • ABYSS -K[K-mer value] input.fastq -o output.fastq

    • perform operation for multiple k-mer value in loop

                 for k in {20..40}; do 

                 ABYSS -k$kreads.fa -o contigs-k$k.fa


    Abyss output file
    ABySS output file

    • ABySS output file consist of contigs generated by the assembly

    • Each contig in the output file consist of 2 lines

    • Line 1 

                >n iii jjj  where n=contig id 

                                       iii=contig length in nucleotides

      jjj=absolute coverage value

    • Line 2

            AAAAACTAATTTCTGAAAT (contig sequence)

    Abyss output
    ABySS output 


    • de Bruijn graph - based assembler

    • best for high coverage very short read data sets

    • leverages paired end information really well


    • velveth - performs hashing

    • default k-mer value - 31  

    • specify input format (-fastq), read-type (-long)


    • velvetg - generates the graph and forms the assembly

    • exp_cov (Expected coverage) - auto  

    Pcap 454

    • Beta version (not yet released)

    • Overlap Layout Consensus assembler

    • Designed for 454 paired-end data

    • Computationally efficient

    Pcap 454 commands
    PCAP-454 commands

    • Input files: Fasta and Qual files separately zipped

    • Fofn file: Specify the name of the input file

    • Run the automated script

    • $./autopcapfofn > auto.log

    Work cited
    Work Cited

    M19107 original reads filtered reads Control with