1 / 81

Outline

Outline. ESP + CGH (continue from last time) CGH analysis using HMM Cancer Progression Models. Array Comparative Genomic Hybridization (aCGH). 4. Combining ESP with other genome data. Joint work with Z. Yakhini, D. Lipson (Agilent and Technion). CGH Analysis.

dunnl
Download Presentation

Outline

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. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Outline • ESP + CGH (continue from last time) • CGH analysis using HMM • Cancer Progression Models

  2. Array Comparative Genomic Hybridization (aCGH) 4. Combining ESP with other genome data Joint work with Z. Yakhini, D. Lipson (Agilent and Technion)

  3. CGH Analysis • Divide genome into segments of equal copy number Copy number profile Genome coordinate Copy number

  4. Copy number Genome coordinate CGH Analysis • Divide genome into segments of equal copy number Copy number profile Numerous methods (e.g. clustering, Hidden Markov Model, Bayesian, etc.) Segmentation • No information about: • Structural rearrangements • (inversions, translocations) • Locations of duplicated material in tumor genome. ESP!

  5. CGH Segmentation How are the copies of segments linked??? 5 3 Copy number 2 Genome Coordinate ES pairs links segments Tumor genome

  6. ESP + CGH ES near segment boundaries 5 3 Copy number 2 Genome Coordinate ESP breakpoint CGH breakpoint

  7. Combining ESP and CGH ES pairs links segments. Copy number balance at each segment boundary: 5 = 2 + 3. 5 3 Copy number 2 Genome Coordinate

  8. Combining ESP and CGH 3 ≤ f(e) ≤ 5 • CGH copy number not exact. • What genome architecture “most consistent” with ESP and CGH data? 5 1 ≤ f(e) ≤ 4 3 1 ≤ f(e) ≤ 3 Copy number 2 Genome Coordinate

  9. Combining ESP and CGH • Edge for each CGH segment. • Edge for each ES pair consistent with segments. • Range of copy number values for each CGH edge. 5 3 Copy number 2 Genome Coordinate 3 ≤ f(e) ≤ 5 1 ≤ f(e) ≤ 3 1 ≤ f(e) ≤ 4 Build graph

  10. Network Flow Problem f(e) Flow constraints: l(e) ≤ f(e) ≤ u(e) CGH edge: l(e) and u(e) from CGH ESP edge: l(e) = 1, u(e) = 1 Flow constraint on each CGH edge l(e) ≤ f(e) ≤ u(e) 8 e

  11. Network Flow Problem f(e) Flow constraints: l(e) ≤ f(e) ≤ u(e) CGH edge: l(e) and u(e) from CGH ESP edge: l(e) = 1, u(e) = 1 Flow in = flow out at each vertex l(e) ≤ f(e) ≤ u(e) 8 e (u,v) f( (u,v) ) = (v,w) f( v,w) ) 8 v

  12. Source/sink Network Flow Problem • Minimum Cost Circulation with Capacity Constraints (Sequencing by Hybridization, Sequence Assembly) Flow constraints: l(e) ≤ f(e) ≤ u(e) CGH edge: l(e) and u(e) from CGH ESP edge: l(e) = 1, u(e) = 1 f(e) Costs: (e) = 0, e ESP or CGH edge 1, e incident to source/sink min e(e) Subject to: l(e) ≤ f(e) ≤ u(e) 8 e (u,v) f( (u,v) ) = (v,w) f( v,w) ) 8 v

  13. Network Flow Results • Unsatisfied flow are putative locations of missing ESP data. • Prioritize further sequencing. f(e) Source/sink • Targeted ESP by screening library with CGH probes.

  14. Network Flow Results • Identify amplified translocations • 14 in MCF7 • 5 in BT474 • Paths of high weight edges: amplicon structures Flow values → Edge weights

  15. Copy number Genome coordinate CGH Analysis • Divide genome into segments of equal copy number Copy number profile Numerous methods (e.g. clustering, Hidden Markov Model, Bayesian, etc.) Segmentation Input: yi = log2 Ti / Ri , clone i = 1, …, N Output: Assignment s(yi) 2 {S1, …, SK} S_k represent copy number states

  16. 1 2 2 1 1 1 1 … 2 2 2 2 … K … … … … x1 K K K K x2 x3 xK … Hidden Markov Models

  17. Example: The Dishonest Casino A casino has two dice: • Fair die P(1) = P(2) = P(3) = P(5) = P(6) = 1/6 • Loaded die P(1) = P(2) = P(3) = P(5) = 1/10 P(6) = 1/2 Casino player switches back-&-forth between fair and loaded die once every 20 turns Game: • You bet $1 • You roll (always with a fair die) • Casino player rolls (maybe with fair die, maybe with loaded die) • Highest number wins $2

  18. Question # 1 – Evaluation GIVEN A sequence of rolls by the casino player 1245526462146146136136661664661636616366163616515615115146123562344 QUESTION How likely is this sequence, given our model of how the casino works? This is the EVALUATION problem in HMMs Prob = 1.3 x 10-35

  19. Question # 2 – Decoding GIVEN A sequence of rolls by the casino player 1245526462146146136136661664661636616366163616515615115146123562344 QUESTION What portion of the sequence was generated with the fair die, and what portion with the loaded die? This is the DECODING question in HMMs FAIR LOADED FAIR

  20. Question # 3 – Learning GIVEN A sequence of rolls by the casino player 1245526462146146136136661664661636616366163616515615115146123562344 QUESTION How “loaded” is the loaded die? How “fair” is the fair die? How often does the casino player change from fair to loaded, and back? This is the LEARNING question in HMMs Prob(6) = 64%

  21. The dishonest casino model 0.05 0.95 0.95 FAIR LOADED P(1|F) = 1/6 P(2|F) = 1/6 P(3|F) = 1/6 P(4|F) = 1/6 P(5|F) = 1/6 P(6|F) = 1/6 P(1|L) = 1/10 P(2|L) = 1/10 P(3|L) = 1/10 P(4|L) = 1/10 P(5|L) = 1/10 P(6|L) = 1/2 0.05

  22. Definition of a hidden Markov model Definition: A hidden Markov model (HMM) • Alphabet = { b1, b2, …, bM } • Set of states Q = { 1, ..., K } • Transition probabilities between any two states aij = transition prob from state i to state j ai1 + … + aiK = 1, for all states i = 1…K • Start probabilities a0i a01 + … + a0K = 1 • Emission probabilities within each state ei(b) = P( xi = b | i = k) ei(b1) + … + ei(bM) = 1, for all states i = 1…K 1 2 K …

  23. A HMM is memory-less At each time step t, the only thing that affects future states is the current state t P(t+1 = k | “whatever happened so far”) = P(t+1 = k | 1, 2, …, t, x1, x2, …, xt) = P(t+1 = k | t) 1 2 K …

  24. 1 1 1 1 … 2 2 2 2 … … … … … K K K K … A parse of a sequence Given a sequence x = x1……xN, A parse of x is a sequence of states  = 1, ……, N 1 2 2 K x1 x2 x3 xK

  25. 1 1 1 1 … 2 2 2 2 … … … … … K K K K … Likelihood of a Parse Simply, multiply all the orange arrows! (transition probs and emission probs) 1 2 2 K x1 x2 x3 xK

  26. 1 1 1 1 … 2 2 2 2 … … … … … K K K K … Likelihood of a parse A compact way to write a01 a12……aN-1N e1(x1)……eN(xN) Number all parameters aij and ei(b); n params Example: a0Fair : 1; a0Loaded : 2; … eLoaded(6) = 18 Then, count in x and  the # of times each parameter j = 1, …, n occurs F(j, x, ) = # parameter j occurs in (x, ) (call F(.,.,.) the feature counts) Then, P(x, ) = j=1…n jF(j, x, ) = = exp[j=1…nlog(j)F(j, x, )] 1 Given a sequence x = x1……xN and a parse  = 1, ……, N, To find how likely is the parse: (given our HMM) P(x, ) = P(x1, …, xN, 1, ……, N) = P(xN, N | N-1) P(xN-1, N-1 | N-2)……P(x2, 2 | 1) P(x1, 1) = P(xN | N)P(N | N-1) ……P(x2 | 2)P(2 | 1)P(x1 | 1)P(1) = a01 a12……aN-1Ne1(x1)……eN(xN) 2 2 K x1 x2 x3 xK

  27. Example: the dishonest casino Let the sequence of rolls be: x = 1, 2, 1, 5, 6, 2, 1, 5, 2, 4 Then, what is the likelihood of  = Fair, Fair, Fair, Fair, Fair, Fair, Fair, Fair, Fair, Fair? (say initial probs a0Fair = ½, aoLoaded = ½) ½  P(1 | Fair) P(Fair | Fair) P(2 | Fair) P(Fair | Fair) … P(4 | Fair) = ½  (1/6)10  (0.95)9 = .00000000521158647211 ~= 0.5  10-9

  28. Example: the dishonest casino So, the likelihood the die is fair in this run is just 0.521  10-9 OK, but what is the likelihood of = Loaded, Loaded, Loaded, Loaded, Loaded, Loaded, Loaded, Loaded, Loaded, Loaded? ½  P(1 | Loaded) P(Loaded, Loaded) … P(4 | Loaded) = ½  (1/10)9  (1/2)1 (0.95)9 = .00000000015756235243 ~= 0.16  10-9 Therefore, it somewhat more likely that all the rolls are done with the fair die, than that they are all done with the loaded die

  29. Example: the dishonest casino Let the sequence of rolls be: x = 1, 6, 6, 5, 6, 2, 6, 6, 3, 6 Now, what is the likelihood  = F, F, …, F? ½  (1/6)10  (0.95)9 = 0.5  10-9, same as before What is the likelihood  = L, L, …, L? ½  (1/10)4  (1/2)6 (0.95)9 = .00000049238235134735 ~= 0.5  10-7 So, it is 100 times more likely the die is loaded

  30. The three main questions on HMMs • Evaluation GIVEN a HMM M, and a sequence x, FIND Prob[ x | M ] • Decoding GIVEN a HMM M, and a sequence x, FIND the sequence  of states that maximizes P[ x,  | M ] • Learning GIVEN a HMM M, with unspecified transition/emission probs., and a sequence x, FIND parameters  = (ei(.), aij) that maximize P[ x |  ]

  31. Let’s not be confused by notation P[ x | M ]: The probability that sequence x was generated by the model The model is: architecture (#states, etc) + parameters  = aij, ei(.) So, P[x | M] is the same with P[ x |  ], and P[ x ], when the architecture, and the parameters, respectively, are implied Similarly, P[ x,  | M ], P[ x,  |  ] and P[ x,  ] are the same when the architecture, and the parameters, are implied In the LEARNING problem we always write P[ x |  ] to emphasize that we are seeking the * that maximizes P[ x |  ]

  32. Problem 1: Decoding Find the most likely parse of a sequence

  33. 1 2 2 1 1 1 1 … 2 2 2 2 … K … … … … x1 K K K K … x2 x3 xK Decoding GIVEN x = x1x2……xN Find  = 1, ……, N, to maximize P[ x,  ] * = argmax P[ x,  ] Maximizes a01 e1(x1) a12……aN-1N eN(xN) Dynamic Programming! Vk(i) = max{1… i-1} P[x1…xi-1, 1, …, i-1, xi, i = k] = Prob. of most likely sequence of states ending at state i = k Given that we end up in state k at step i, maximize product to the left and right

  34. Decoding – main idea Inductive assumption: Given that for all states k, and for a fixed position i, Vk(i) = max{1… i-1} P[x1…xi-1, 1, …, i-1, xi, i = k] What is Vl(i+1)? From definition, Vl(i+1) = max{1… i}P[ x1…xi, 1, …, i, xi+1, i+1 = l ] = max{1… i}P(xi+1, i+1 = l | x1…xi, 1,…, i) P[x1…xi, 1,…, i] = max{1… i}P(xi+1, i+1 = l | i ) P[x1…xi-1, 1, …, i-1, xi, i] = maxk [P(xi+1, i+1 = l | i=k) max{1… i-1}P[x1…xi-1,1,…,i-1, xi,i=k]] = maxk [ P(xi+1 | i+1 = l ) P(i+1 = l | i=k) Vk(i)] = el(xi+1)maxk aklVk(i)

  35. The Viterbi Algorithm Input: x = x1……xN Initialization: V0(0) = 1 (0 is the imaginary first position) Vk(0) = 0, for all k > 0 Iteration: Vj(i) = ej(xi)  maxk akj Vk(i – 1) Ptrj(i) = argmaxk akj Vk(i – 1) Termination: P(x, *) = maxk Vk(N) Traceback: N* = argmaxk Vk(N) i-1* = Ptri (i)

  36. The Viterbi Algorithm x1 x2 x3 ………………………………………..xN Similar to “aligning” a set of states to a sequence Time: O(K2N) Space: O(KN) State 1 2 Vj(i) K

  37. Viterbi Algorithm – a practical detail Underflows are a significant problem P[ x1,…., xi, 1, …, i ] = a01 a12……ai e1(x1)……ei(xi) These numbers become extremely small – underflow Solution: Take the logs of all values Vl(i) = logek(xi) + maxk [ Vk(i-1) + log akl ]

  38. Example Let x be a long sequence with a portion of ~ 1/6 6’s, followed by a portion of ~ ½ 6’s… x = 123456123456…123456626364656…1626364656 Then, it is not hard to show that optimal parse is (exercise): FFF…………………...FLLL………………………...L 6 characters “123456” parsed as F, contribute .956(1/6)6 = 1.610-5 parsed as L, contribute .956(1/2)1(1/10)5 = 0.410-5 “162636” parsed as F, contribute .956(1/6)6 = 1.610-5 parsed as L, contribute .956(1/2)3(1/10)3 = 9.010-5

  39. Problem 2: Evaluation Find the likelihood a sequence is generated by the model

  40. 1 1 1 1 … 2 2 2 2 … … … … … K K K K … Generating a sequence by the model Given a HMM, we can generate a sequence of length n as follows: • Start at state 1 according to prob a01 • Emit letter x1 according to prob e1(x1) • Go to state 2 according to prob a12 • … until emitting xn 1 a02 2 2 0 K e2(x1) x1 x2 x3 xn

  41. A couple of questions P(box: FFFFFFFFFFF) = (1/6)11 * 0.9512 = 2.76-9 * 0.54 = 1.49-9 P(box: LLLLLLLLLLL) = [ (1/2)6 * (1/10)5 ] * 0.9510 * 0.052 = 1.56*10-7 * 1.5-3= 0.23-9 Given a sequence x, • What is the probability that x was generated by the model? • Given a position i, what is the most likely state that emitted xi? Example: the dishonest casino Say x = 12341…23162616364616234112…21341 Most likely path:  = FF……F However: marked letters more likely to be L than unmarked letters F F

  42. Evaluation We will develop algorithms that allow us to compute: P(x) Probability of x given the model P(xi…xj) Probability of a substring of x given the model P(I = k | x) Probability that the ith state is k, given x A more refined measure of which states x may be in

  43. The Forward Algorithm We want to calculate P(x) = probability of x, given the HMM Sum over all possible ways of generating x: P(x) =  P(x, ) =  P(x | ) P() To avoid summing over an exponential number of paths , define fk(i) = P(x1…xi, i = k) (the forward probability)

  44. The Forward Algorithm – derivation Define the forward probability: fk(i) = P(x1…xi, i = k) = 1…i-1 P(x1…xi-1, 1,…, i-1, i = k) ek(xi) = l1…i-2 P(x1…xi-1, 1,…, i-2, i-1 = l) alk ek(xi) = lP(x1…xi-1, i-1 = l) alk ek(xi) = ek(xi) lfl(i-1) alk

  45. The Forward Algorithm We can compute fk(i) for all k, i, using dynamic programming! Initialization: f0(0) = 1 fk(0) = 0, for all k > 0 Iteration: fk(i) = ek(xi) l fl(i-1) alk Termination: P(x) = k fk(N) ak0 Where, ak0 is the probability that the terminating state is k (usually = a0k)

  46. Relation between Forward and Viterbi FORWARD Initialization: f0(0) = 1 fk(0) = 0, for all k > 0 Iteration: fl(i) = el(xi) k fk(i-1) akl Termination: P(x) = k fk(N) ak0 VITERBI Initialization: V0(0) = 1 Vk(0) = 0, for all k > 0 Iteration: Vj(i) = ej(xi) maxkVk(i-1) akj Termination: P(x, *) = maxkVk(N)

  47. Motivation for the Backward Algorithm We want to compute P(i = k | x), the probability distribution on the ith position, given x We start by computing P(i = k, x) = P(x1…xi, i = k, xi+1…xN) = P(x1…xi, i = k) P(xi+1…xN | x1…xi, i = k) = P(x1…xi, i = k) P(xi+1…xN | i = k) Then, P(i = k | x) = P(i = k, x) / P(x) Forward, fk(i) Backward, bk(i)

  48. The Backward Algorithm – derivation Define the backward probability: bk(i) = P(xi+1…xN | i = k) = i+1…N P(xi+1,xi+2, …, xN, i+1, …, N | i = k) = li+1…N P(xi+1,xi+2, …, xN, i+1 = l, i+2, …, N | i = k) = l el(xi+1) akli+1…N P(xi+2, …, xN, i+2, …, N | i+1 = l) = l el(xi+1) aklbl(i+1)

  49. The Backward Algorithm We can compute bk(i) for all k, i, using dynamic programming Initialization: bk(N) = ak0, for all k Iteration: bk(i) = l el(xi+1) akl bl(i+1) Termination: P(x) = l a0l el(x1) bl(1)

  50. Computational Complexity What is the running time, and space required, for Forward, and Backward? Time: O(K2N) Space: O(KN) Useful implementation technique to avoid underflows Viterbi: sum of logs Forward/Backward: rescaling at each position by multiplying by a constant

More Related