Sequence Assembly and Protein Docking Algorithms

1. Sequence Assembly and Protein Docking Algorithms Vicky Choi Department of Computer Science Duke University

2. Outline • Sequence Assembly Algorithm (Joint work with Martin Farach-Colton @Rutgers) • Local Search Algorithm for Rigid Protein Docking (Joint work with Pankaj K. Agarwal, Herbert Edelsbrunner, Johannes Rudolph @Duke)

3. Outline: Sequence Assembly • Biological Background • Human Genome Project and the Sequence Assembly Problem • The BARNACLEAssembler

4. DNA A DNA molecule consists of two strands which are tied together in a helical structure. Image Credit: US Department of Energy Human Genome Program http://www.ornl.gov/hgmis Each strand is represented by a string over the alphabet {A,C,G,T}, called a DNA sequence. Example AAGCTTCAGTTCCTGACCTTCCAATCGCAA {A,C,G,T} = nucleotide, base, basepair (bp)

5. (3’) (5’) one strand ) another strand Example (5') (3') ACCATGGTGCACCTGACTCCTGAGGAG TGGTACCACGTGGACTGAGGACTCCTC (3') (5') (3’) (5’) Two Strands: Reverse Complementary Orientation: 5’! 3’ Complement: A $ T, C $ G Image Credit: US Department of Energy Human Genome Program http://www.ornl.gov/hgmis

6. A genome is the complete set of DNA sequences of an organism. Human Chromosomes Image Credit: Sanger Center http://www.sanger.ac.uk/ Human Genome ~ 3x109 bp

7. C G A A T C G T C G A T G C T A A T G DNA Sequencing DNA Sequencing is the process for determining the sequence of nucleotides of a region of DNA. Current technology : ~500bp Question: How to sequence a longer stretch of DNA?

8. DNA Cloning Copies of Target Shotgun ACGTAAGAGTACCGATTGGCCA DNA Sequencing Sequence Reads Assembly Consensus Directed Read Final …ACGTAGTCTTAGATGATAGTAGA… Shotgun Sequencing Target DNA

9. Shotgun Sequencing History • 1980s: 5 to 10 Kbp • 1990: 40 Kbp • 1995: 1.8 Mbp (H. Influenzae) • 2000: draft Drosophila (120 Mbp) • 2001: draft Human Genome (3x109bp) (attempted by Celera)

10. Outline: Sequence Assembly • Biological Background • Human Genome Project and the Sequence Assembly Problem • The BARNACLE Assembler

11. Human Genome Project (HGP) • 1988: “Mapping and Sequencing the Human Genome” • 1990: HGP started in US • 2001: A “working draft” version • 2003: Completed by HGP Consortium standard

12. Physical Map Tiling Path of BACs ShotgunSequence & Assemble of each TP BAC Final Sequence Hierarchical Shotgun Sequencing(BAC-by-BAC) • Map First, Then Sequence Human Genome A BAC is a segment of DNA from a chromosome. Each BAC is ~100-200Kb. BAC library (100-200Kb)

13. Human Genome BAC library (100-200Kb) Physical Map ? ? ? ? Physical Map is difficult to build! (original expected time: 5 years) Hierarchical Shotgun Sequencing(BAC-by-BAC) • Map First??

14. BAC Sequence Reads Fragments Phase 1: Draft Ordered Fragments Phase 2: Draft Phase 3: Finished BAC-by-BAC! BAC-Based New Idea: Map + Sequencing concurrently Randomly pick BACs (not wait for Physical Map) and shotgun sequence BACs

15. BAC-Based Human Genome BAC library (100-200Kb) Finished + Draft BACs Sequence Assembly Problem the working draft of the human genome

16. Outline: Sequence Assembly • Biological Background • Human Genome Project and the Sequence Assembly Problem • The BARNACLE Assembler • Details of Input • Difficulties • Basic Idea

17. Details of Input • Sequence Information: • BACs • Overlap Information: • Local Alignments • NT-pairs • Orientation Information: Plasmid, EST, mRNA

18. BAC fragment Input: Sequence Information Recall: A BAC is a contiguous stretch of DNA from a chromosome. Each comes as a set of fragments. • Phase 1,2 = Draft • Phase 3 = Finished

19. Input: Overlap Information • Preprocessing: • Local alignments of all fragment pairs • NT-pairs: Generated from GenBank annotation • submitted from genome centers

20. Example: Input of Dec 2001 freeze Sequence Information: Chromosome Assignments: 31543 by STS; 2450 by Genbank; 1083 unknown Overlap Information: 403,466 fragment pairs, 12,656 NT-pairs Orientation Information: 321,751 fragment pairs

21. True Overlap Repeat-induced Overlap True vs Repeat-induced Overlap

22. Repeats of the Human Genome • High-copy repeats e.g. ALU, L1 • Low-copy repeats (segmental duplication) • Large block (>200Kb) • Highly Similar (>97%)

23. Noise • False positives (FP): due to repeat • False negatives (FN): polymorphism, draft quality • Chimeric BAC (CB)

24. The Basic Idea 1. “Conservatively” assemble fragments

25. A overlaps with B B B B A A A A A overlaps with C C C C Yes. No. Necessary Condition for True Overlaps Does B overlap with C? Idea: assemble non-conflict overlaps first

26. BAC Graph The Basic Idea 1. “Conservatively” assemble fragments into subcontigs

27. Interval Graph The BAC graph is an interval graph! Definition: A graph G is called an interval graph if there is one-one correspondence between its vertices and a set of intervals on the real line such that two vertices are adjacent in G iff their corresponding intervals overlap.

28. Under-represented Necessary… But Not Sufficient Long Repeats

29. Non-interval Graph Collapsing Repeats: Chimeric BAC

30. Forbidden Subgraphs Theorem (Lekkerkerker & Boland 1962) A graph is interval iff it does not contain one of the following (induced) subgraph:

31. Resolving Non-interval Graphs Definition: A vertex u 2 V is I-critical if G|V\{u} is interval. Given a non-interval graph G, identify a forbidden subgraph. If at least one of the vertices of the forbidden subgraph is I-critical, we say G is fixable. • Based on the structure of the forbidden subgraph, • a fixable graph G is resolved by • adding an FN edge; or • removing FP edges; or • removing a vertex.

32. Divide and Conquer Method For the non-fixable graphs, we employ a divide-and-conquer method by dividing the graph according to some articulation points such that each subcomponent is fixable.

33. BAC Graph The Basic Idea 1. “Conservatively” assemble fragments into subcontigs 2. Resolve Non-interval Graph and Find an Interval Realization of the BAC Graph 3. Orient and order subcontigs

34. Error Detection 1. “Conservatively” assemble fragments into subcontigs • wrong NT-pairs (annotation from genome centers) • chromosome misassignments 2. Resolve Non-interval Graph and Find an Interval Realization of the BAC Graph • chimerics 3. Orient and order subcontigs • fragment misassemblies This is the only algorithm available that does Error Detection.

35. Output: Contigs

36. Other Two Assemblers • GigAssembler by Jim Kent and David Haussler • (stop after April 2001 freeze) • NCBI’s assembler – top-down approach: • build a physical map using sequence overlaps as • fingerprint overlaps; • using some scoring functions to resolve conflicts.

37. NCBI’s assembly BARNACLE’s assembly

38. Comparison with NCBI’s Assembly(Dec 2001)

39. How was the human genome “finished”? • Hand-curate tiling path of BACs (by Genome Centers) • Finish sequencing the tiling path of BACs only • Assemble by NCBI’s assembler based on the • hand-curated tiling paths

40. Incorporating segmental duplication database Collaboration with Evan E. Eichler (Department of Genetics, Case Western Reserve University) • BARNACLE’s assembly suggested that at least 89 repeat-contained BACs were dropped from the tiling path. • 69 were added to HGP’s final tiling path • 20 were declared unnecessary • Due to disagreement about repeat structure of genome The Sequence and Assembly of Highly Duplicated Regions in the Human Genome. V. Choi, J. Bailey, G. Schuler, Z. Gu, P. Li, M. Farach-Colton and E. Eichler. Genome Sequencing & Biology meeting at Cold Spring Harbor Laboratory 2002.

41. Conclusions • Better assembly • Error detection • Measured by the assembled BAC length • Efficient (3 minutes on a Pentium III) Reference : V. Choi, M. Farach-Colton. BARNACLE: An Assemble Algorithm for Clone-based Sequences of Whole Genomes. Gene, 320, 165-176, 2003. • To do large scale sequencing: • Handle repeats • Design in data acquisition that will permit error detection & correction

42. Acknowledgement • Wojciech Makalowski (NCBI/Penn State University) • David Lipman (NCBI) • Greg Schuler (NCBI) • Evan E. Eichler (Case Western Reserve University) • Granger Sutton (Celera) • JinSheng Lai (Waksman Institute, Rutgers University) • NCBI/NIH pre-doctoral visiting fellowship • Program in Mathematics and Molecular Biology (PMMB), • Burroughs Wellcome Fund Interface Program fellowship

43. Outline : Protein Docking • Protein-Protein Docking • Local Search Algorithm • Test Results

44. Protein-Protein Docking Barstar Barnase 1BRS : Barnase + Barstar

45. Protein Re-Docking Problem(Bound Protein Docking) Given a known protein-protein complex A-B(native configuration), randomly separate two proteins. Fix A, find a rigid motion m such that m(B) is near-native. Rigid Body Assumption

46. Formulation of Rigid Protein Docking • A scoring function that can discriminate correct docking configuration from incorrect ones; • A search algorithm that finds the docking configuration measured by the scoring function.

47. Protein A protein molecule consists of a set of atoms. Each atom is represented by a ball in R3. Notation: A = { (a1, r1), …, (an,rn) } where ai2R3 is the ith atom center with (van der Waals) radius ri

48. Our Scoring Function

49. Exhaustive Search • Sampling the rigid motion space (6-dimension) • Evaluate each motion using the scoring function Rigid Motion = Rotation + Translation A rotation in R3 can be specified by a rotation angle q about a rotation axis u – represented by unit quaternion. Sampling Rotation Space ) S3 (unit sphere in R4) A translation in R3 is a 3-dimensional vector (x,y,z) 2R3. Sampling Translation Space : a 3-dimensional grid

50. Protein Re-Docking Without False Positives Empirical Results: • The configuration m(B) for which • Score(A,m(B)) is maximized; • Bump(A,m(B)) · 7; • isnear-native :RMSD(Bnative, m(B)) · 3. Other prior works (e.g. FFT-based, geometric hashing) generate multiple possible docking configurations (i.e. near-native + false positives).