MW 11:00-12:15 in Beckman B302 Prof: Gill Bejerano TAs: Jim Notwell & Harendra Guturu

1. CS173 Lecture 12: Chains & Nets, Conservation & Function MW  11:00-12:15 in Beckman B302 Prof: Gill Bejerano TAs: Jim Notwell & Harendra Guturu http://cs173.stanford.edu [BejeranoWinter12/13]

2. Announcements • HW2 Due Today • As are project assignments • Coming monday 2/25 lecture has been moved to LK101(building next door – we’ll post instructions) http://cs173.stanford.edu [BejeranoWinter12/13]

3. Inferring Genomic Mutations From Alignments of Genomes http://cs173.stanford.edu [BejeranoWinter12/13]

4. Gene tree Speciation Duplication Loss Terminology Orthologs : Genes related via speciation (e.g. C,M,H3) Paralogs: Genes related through duplication (e.g. H1,H2,H3) Homologs: Genes that share a common origin (e.g. C,M,H1,H2,H3) single ancestral gene Species tree http://cs173.stanford.edu [BejeranoWinter12/13]

5. What? • Compare whole genomes • Compare two genomes • Within (intra) species • Between (inter) species • Compare genome to itself • Why? • Comparison reveals functional and neutral regions • Homologous regions most often have similar functions • Modification of functional regions can reveal • Disease susceptibility • Adaptation • How? http://cs173.stanford.edu [BejeranoWinter12/13]

6. Every Genome is Different DNA Replication is imperfect – between individuals of the same species, even between the cells of an individual. junk functional ...ACGTACGACTGACTAGCATCGACTACGA... chicken TT CAT egg ...ACGTACGACTGACTAGCATCGACTACGA... many changes are not tolerated “anything goes” chicken http://cs173.stanford.edu [BejeranoWinter12/13]

7. Sequence Alignment AGGCTATCACCTGACCTCCAGGCCGATGCCC TAGCTATCACGACCGCGGTCGATTTGCCCGAC -AGGCTATCACCTGACCTCCAGGCCGA--TGCCC--- TAG-CTATCAC--GACCGC--GGTCGATTTGCCCGAC Similarity is often measured using “%id”, or percent identity %id = number of matching bases / number of alignment columns Where Every alignment column is a match / mismatch / indel base Where indel = insertion or deletion (requires an outgroup to resolve)

8. What to expect from genome comparisons? human Objective: find local alignment blocks, that are likely homologous (share common origin) O(mn) examine the full matrix using DP O(m+n) heuristics based on seeding + extension trades sensitivity for speed lizard http://cs173.stanford.edu [BejeranoWinter12/13]

9. “Raw” (B)lastz track (no longer displayed) Alignment = homologous regions Protease Regulatory Subunit 3

10. Chaining co-linear alignment blocks human Objective: find local alignment blocks, that are likely homologous (share common origin) Chaining strings together co-linear blocks in the target genome to which we are comparing. Double lines when there is unalignable sequence in the other species. Single lines when there isn’t. lizard http://cs173.stanford.edu [BejeranoWinter12/13]

11. Reference genome perspective,The Use of an Outgroup Outgroup Sequence A B C Human Sequence Mouse Sequence D E A B C A B C D E B’ D E In Human Browser In Mouse Browser Implicit Human sequence Implicit Mouse sequence … … D E … … Mouse chains Human chains D E D E B’

12. Gap Types: Single vs Double sided Ancestral Sequence A B C D E Human Sequence Mouse Sequence A B C A B C D E B’ D E In Human Browser In Mouse Browser Implicit Human sequence Implicit Mouse sequence … … D E … … Mouse chains Human chains D E D E B’

13. Conservation Track Documentation http://cs173.stanford.edu [BejeranoWinter12/13]

14. Chains join together related local alignments likely ortholog likely paralogs shared domain? Protease Regulatory Subunit 3 http://cs173.stanford.edu [BejeranoWinter12/13]

15. Note: repeats are a nuisance human If, for example, human and mouse have each 10,000 copiesof the same repeat: We will obtain and need to output 108 alignments of all these copies to each other. Note that for the sake of this comparison interspersed repeats and simple repeats are equal nuisances. Also note that simple repeats, but not interspersed repeats, violate the assumption that similar sequences are homologous. mouse Solution: 1 Discover all repetitive sequences in each genome. 2 Mask them when doing genome to genome comparison. 3 Chain your alignments. 4 Add back to the alignments only repeat matches that lie within pre-computed chains. http://cs173.stanford.edu [BejeranoWinter12/13]

16. Chains • a chain is a sequence of gapless aligned blocks, where there must be no overlaps of blocks' target or query coords within the chain. • Within a chain, target and query coords are monotonically non-decreasing. (i.e. always increasing or flat) • double-sided gaps are a new capability (blastz can't do that) that allow extremely long chains to be constructed. • not just orthologs, but paralogs too, can result in good chains. but that's useful! • chains should be symmetrical -- e.g. swap human-mouse -> mouse-human chains, and you should get approx. the same chains as if you chain swapped mouse-human blastz alignments. • chained blastz alignments are not single-coverage in either target or query unless some subsequent filtering (like netting) is done. • chain tracks can contain massive pileups when a piece of the target aligns well to many places in the query. Common causes of this include insufficient masking of repeats and high-copy-number genes (or paralogs). [Angie Hinrichs, UCSC wiki] http://cs173.stanford.edu [BejeranoWinter12/13]

17. Before and After Chaining http://cs173.stanford.edu [BejeranoWinter12/13]

18. Chaining Algorithm Input - blocks of gapless alignments from blastz Dynamic program based on the recurrence relationship:score(Bi) = max(score(Bj) + match(Bi) - gap(Bi, Bj)) Uses Miller’s KD-tree algorithm to minimize which parts of dynamic programming graph to traverse. Timing is O(N logN), where N is number of blocks (which is in hundreds of thousands) j<i See [Kent et al, 2003] “Evolution's cauldron: Duplication, deletion, and rearrangement in the mouse and human genomes” http://cs173.stanford.edu [BejeranoWinter12/13]

19. Netting Alignments Commonly multiple mouse alignments can be found for a particular human region, particularly for coding regions. Net finds best match mouse match for each human region. Highest scoring chains are used first. Lower scoring chains fill in gaps within chains inducing a natural hierarchy. http://cs173.stanford.edu [BejeranoWinter12/13]

20. Net highlights rearrangements A large gap in the top level of the net is filled by an inversion containing two genes. Numerous smaller gaps are filled in by local duplications and processed pseudo-genes. http://cs173.stanford.edu [BejeranoWinter12/13]

21. A Rearrangement Hot Spot Rearrangements are not evenly distributed. Roughly 5% of the genome is in hot spots of rearrangements such as this one. This 350,000 base region is between two very long chains on chromosome 7. http://cs173.stanford.edu [BejeranoWinter12/13]

22. Nets attempt to capture the ortholog (they also hide everything else) http://cs173.stanford.edu [BejeranoWinter12/13]

23. Retroposed Genes and Pseudogenes Pseudogenes (“dead genes”): Genomic sequences that resemble (originated from) genes that no longer make proteins. Retrogenes (“retrotranscribed”): Protein coding RNA that was reverse transcribed and inserted back into the genome. The RNA can be grabbed at any stage (partial/full transcript, before/during/after all introns are spliced). http://cs173.stanford.edu [BejeranoWinter12/13]

24. Useful in finding pseudogenes gene pred. Ensembl and Fgenesh++ automatic gene predictions confounded by numerous processed pseudogenes. Domain structure of resulting predicted protein must be interesting! http://cs173.stanford.edu [BejeranoWinter12/13]

25. Nets/chains can reveal retrogenes (and when they jumped in!) http://cs173.stanford.edu [BejeranoWinter12/13]

26. Nets • a net is a hierarchical collection of chains, with the highest-scoring non-overlapping chains on top, and their gaps filled in where possible by lower-scoring chains, for several levels. • a net is single-coverage for target but not for query. • because it's single-coverage in the target, it's no longer symmetrical. • the netter has two outputs, one of which we usually ignore: the target-centric net in query coordinates. The reciprocal best process uses that output: the query-referenced (but target-centric / target single-cov) net is turned back into component chains, and then those are netted to get single coverage in the query too; the two outputs of that netting are reciprocal-best in query and target coords. Reciprocal-best nets are symmetrical again. • nets do a good job of filtering out massive pileups by collapsing them down to (usually) a single level. • GB: for human inspection always prefer looking at the chains! [Angie Hinrichs, UCSC wiki] http://cs173.stanford.edu [BejeranoWinter12/13]

27. Before and After Netting http://cs173.stanford.edu [BejeranoWinter12/13]

28. Convert / LiftOver "LiftOver chains" are actually chains extracted from nets, or chains filtered by the netting process. LiftOver – batch utility http://cs173.stanford.edu [BejeranoWinter12/13]

29. What nets can’t show, but chains will http://cs173.stanford.edu [BejeranoWinter12/13]

30. Same Region… same in all the other fish http://cs173.stanford.edu [BejeranoWinter12/13]

31. Drawbacks • Inversions not handled optimally Chains > > > > chr1 > > > > > > > chr1 > > > < < < < chr5 < < < < < < < < chr1 < < < < Nets > > > > chr1 > > > > > > > chr1 > > > < < < < chr5 < < < <

32. Drawbacks • High copy number genes can break orthology

33. Self Chain reveals paralogs (self net is meaningless) http://cs173.stanford.edu [BejeranoWinter12/13]

34. Conservation and Function http://cs173.stanford.edu [BejeranoWinter12/13]

35. Evolution = Mutation + Selection Mistakes can happen during DNA replication. Mistakes are oblivious to DNA segment function. But then selection kicks in. junk functional ...ACGTACGACTGACTAGCATCGACTACGA... chicken TT CAT egg ...ACGTACGACTGACTAGCATCGACTACGA... many changes are not tolerated “anything goes” chicken Conservation implies function! (But what function?) http://cs173.stanford.edu [BejeranoWinter12/13]

36. Vertebrates: what to sequence? , Stickleback Which species to compare to? Too close and purifying selection will be largely indistinguishable from the neutral rate. Too far and many functional orthologs will diverge beyond our ability to accurately align them. , Lizard , Opossum too far sweet spot too close  you are here [Human Molecular Genetics, 3rd Edition] http://cs173.stanford.edu [BejeranoWinter12/13]

37. Searching Near And Far Search too near (eg human to chimp or orang above) and you cannot distinguish neutral sequence from sequence under purifying selection. Search further still (eg mouse) and the two distributions pry apart.But now you’ve lost younger functional sequences born after the split. Ie, conservation implies function, butlack of conservation does NOT imply lack of function! http://cs173.stanford.edu [BejeranoWinter12/13]

38. Phylogenetic Shadowing , Stickleback , Lizard , Opossum “too close” can actually be a boon if you have enough closely related genomes too close  you are here http://cs173.stanford.edu [BejeranoWinter12/13]

39. PhastCons Conserved Elements http://cs173.stanford.edu [BejeranoWinter12/13]

40. Distant homologies When species diverge too much (e.g. chicken and beyond above), confident alignments can no longer be detected at the DNA level. E.g.: all SPI1 and SLC39A13 exons are there in chicken & fish. http://cs173.stanford.edu [BejeranoWinter12/13]

41. Distant homologies search strategies Here it is much better to search a gene model from species A (e.g human) against the genome of species B (e.g. chicken) This is a search of amino acids in all their possible codons into a gene structure with unknown exon – intron structure. (eg TBLASTN, translated BLAT) http://cs173.stanford.edu [BejeranoWinter12/13]

42. Distant homologies Find the most distantly related genes using gene models in both species: 1 search amino acids sequences against each other. (eg using BLASTP). 2 Map your hits back to the two respective genomes, anchored on the amino acid alignment (respecting any exon-intron gene body structure change). 3 Examine co-linear homology of flanking genes to try and call orthologs from paralogs. http://cs173.stanford.edu [BejeranoWinter12/13]

43. RNA homology searches 1 Define a mathematical construct that describe potential homologs. 2 Go search for them (efficiently!). 3 Examine genomic context. http://cs173.stanford.edu [BejeranoWinter12/13]

44. Enhancer remote homologs Enhancer = • Gene regulatory sequences in general are the most challenging to search for: • Individual binding sites are very flexible. • Gaps between binding sites may evolve (semi) neutrally,making DNA alignment seeding particularly frail. • Binding site gain/loss and shuffling may or may not be allowed – we need a better understanding of underlying logic. http://cs173.stanford.edu [BejeranoWinter12/13]

45. Exceptionally Old Enhancers Exist But how many of these really exist? http://cs173.stanford.edu [BejeranoWinter12/13]

46. * * * * * Ultraconservation: No known function requires this much conservation ? CDS ncRNA TFBS seq. http://cs173.stanford.edu [BejeranoWinter12/13]

47. “Gene” Finding III: Comparative Genomics http://cs173.stanford.edu [BejeranoWinter12/13]

48. The challenge: map code to output 1013 cells 3*109 letters genome person • Ultimately we sequence genomes, and study their function in detail to understand genome to phenotype relationships: • Minus side: Genomic contribution to disease • Plus side: Adaptation and speciation To be continued… http://cs173.stanford.edu [BejeranoWinter12/13]