Luận văn Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - áp dụng trên cây sinh loài nhỏ

ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 27 Công thức này thoạt nhìn không có vẻ đơn giản, nhưng thật ra nó rất dễ dàng khi tính đạo hàm cho log likelihood trong trường hợp này cũng như nhiều trường hợp khác. Lấy đạo hàm và cho chúng bằng 0, chúng ta được: ' (1 ) ( )( ) 0 1 .(1 ) .(1 ) H T H T H H T D N N N N N N Nl θ θ θθ θ θ θ θ θ θ − − − += − = = =− − − TH H NN N +=⇔θ  Bảng 1: Bảng biến thiên của hàm hợp lý với θ là nghiệm chúng ta cần tìm, phù hợp với những gì chúng ta mong muốn. Theo ví dụ trên nếu (NH, NT ) = (3, 2) và MLE tính được là 3 0.6 5 = . Đồ thị của hàm hợp lý cho ta thấy ở hình 2. Hình 2: Đồ thị của hàm hợp lý 3.3. Tổng quát hóa bài toán ước lượng hợp lý cực đại 3.3.1. Ước lượng hợp lý cực đại trên mẫu quan sát Nếu x là biến ngẫu nhiên với hàm phân bố: [ ] 1 2( , , ..., )x i Kf θ θ θ với 1 2, , ..., Kθ θ θ là K tham số cần phải ước lượng, với dãy N mẫu độc lập là x[1], x[2], ..., x[N]. Thì hàm likelihood được cho bởi tích sau: 1 2 [ ] 1 2 1 ( , , ..., ) ( , , ..., ) N D K x i i L f Kθ θ θ θ θ θ = =∏ Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 28 và hàm ln likelihood như sau: [ ] 1 1 1 ( ) ln ( ) ln ( , , ..., ) N D D x i k i L fθ θ θ θ = = =∑A θ MLE của 1 2, , ..., Kθ θ θ đạt được khi ( )DL θ hay ( )D θA là lớn nhất, chúng ta đã biết xác định giá trị lớn nhất với ( )D θA dễ hơn với ( )DL θ , vậy MLE của 1 2, , ..., Kθ θ θ là giải hệ K phương trình sau: ( ) 0, 1,2, ..., j j Kθ ∂ = =∂ A Ví dụ: Tung một con xúc sắc có K = 6 mặt, chúng ta muốn xác định những tham số 1 2, ,..., Kθ θ θ là xác suất của mặt có nút tương ứng 1, 2,…, K nhận được khi tung xúc sắc. Từ quan sát ta có 1 2, ,..., KN N N là số lượng tương ứng của từng mặt khi quan sát. Theo công thức hàm khả năng sẽ: 1 ( ) j K N D j j L θ θ = =∏ và hàm ln likelihood tương ứng sẽ là: 1 ( ) ln ( ) ln( ) K D D i i l L N iθ θ θ = = =∑ Sau khi giải hệ phương trình ( ) 0, 1,2,..., j l j Kθ ∂ = =∂ chúng ta được: 1 k k K l l N N θ = = ∑  với k =1, …, K 3.3.2. Một số phương pháp giải phương trình hợp lý Theo trên, giải phương trình hợp lý làm cực đại phương trình: 1 2( ) ( ) ( ) ... ( )1 2 u u nL f f f n uθ θ θ θ= với iu ∈` . Hiện nay có hai hướng tiếp cận khác nhau để giải quyết bài toán này, trong mỗi phương pháp có những ưu và khuyết điểm riêng của nó: Phương pháp gần đúng: Giải phương trình hợp lý bằng phương pháp tìm kiếm cục bộ, heuristics, …Ưu điểm của phương pháp này là nhanh chóng, có thể giải quyết Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 29 trên những bài toán lớn. Nhược điểm lớn nhất của phương pháp này là tính tin cậy không cao. Phương pháp tính toán đại số: Ngược lại với phương pháp gần đúng trên, phương pháp tính toán đại số hiện nay chỉ giải quyết được với những bài toán nhỏ, nhưng cho kết quả chính xác. Với sự tiến bộ của khoa học kỹ thuật nói chung và ngành máy tính cũng như lãnh vực đại số máy tính nói riêng, đã mở ra con đường cho hướng tiếp cận này. Vì lý do trên phương pháp này được chọn sử dụng để giải quyết bài toán ước lượng hợp lý cực đại - áp dụng trên cây sinh loài nhỏ. Để hiểu rõ cây sinh loài, ước lượng hợp lý cực đại trên cây sinh loài chúng ta tìm hiểu sơ qua cây sinh loài và mô hình xác suất thống kê trên cây sinh loài ở chương sau. Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 30 Chương 4. CÂY SINH LOÀI - MÔ HÌNH XÁC SUẤT THỐNG KÊ TRÊN CÂY SINH LOÀI Chương này giới thiệu cây sinh loài cũng mô hình xác suất thống kê trên cây sinh loài. Ngoài ra cũng giới thiệu một số mô hình thường sử dụng hiện nay trên cây sinh loài như mô hình Neyman 2 trạng thái, Jukes – Cantor, Kimura với 2 và 3 tham số. 4.1. Giới thiệu sơ lược về cây sinh loài Cây sinh loài (còn gọi là cây tiến hóa hay là cây chủng loài) mô tả lịch sử tiến hóa của một nhóm các loài (species) với những đặc tính khác nhau nhưng cùng có mối quan hệ họ hàng với nhau và cùng hình thành từ một tổ tiên chung trong quá khứ. Có nhiều hướng nghiên cứu khác nhau để chứng minh đặc điểm phát sinh sinh loài này. Trước hết, người ta có thể so sánh trình tự các đoạn DNA (thuộc sinh học phân tử hay hệ gene học (genomics); hoặc so sánh các hóa thạch (fossil) hoặc các di chỉ (record) của sinh vật cổ (thuộc khảo cổ học - paleontology). Các nhà sinh học tổ chức và phân tích các mối quan hệ tiến hóa thông qua các phương pháp khác nhau, bao gồm phân loại học (phylogenetics), ngoại hình học (phenetics) và cladistics. Các sự kiện chính xảy ra trong quá trình tiến hóa của sự sống được xây dựng thành biểu đồ thời gian của tiến hóa (evolutionary timeline) dựa trên các hiểu biết hiện nay của khoa học. Hình 3 cho ta thấy hình dạng của cây sinh loài sự sống trên hành tinh chúng ta. Hình 3: Cây sinh loài của sự sống Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 31 4.2. Các nghiên cứu phát sinh sinh loài Trong ngành sinh học, người ta nghiên cứu mối quan hệ giữa các loài sinh vật thông qua các bằng chứng phân tử, cụ thể là trình tự DNA và protein. Như vậy sự khác biệt giữa các trình tự (DNA) chỉ định sự phân kỳ di truyền như là kết quả của tiến hóa phân tử theo tiến trình thời gian. Các phương pháp dùng để nghiên cứu phát sinh sinh loài chủ yếu dựa trên một sự giả định về các tiến trình tiến hóa ở mức phân tử thông qua việc quan sát phân tích trình tự DNA hoặc protein. Bằng cách sử dụng công cụ máy tính, các chuỗi dữ liệu sẽ được mô phỏng tiến trình tiến hóa và phân tích tiến trình phát sinh sinh loài. Giả sử là chúng ta có một “cây tiến hóa đúng”, chúng ta có thể dùng nó để kiểm tra lại độ chính xác, tính nhất quán khả năng tin cậy của những mô hình tiến hóa. Tuy nhiên khi sử dụng các dữ liệu sinh học, cái gọi là cây tiến hóa có thể không bao giờ có, hoặc ít ra cũng có thể nói là KHÔNG BIẾT. Do vậy người ta chấp nhận một cây tiến hóa được dựng nên mà người ta tin là nó GIỐNG NHẤT với cây tiến hóa đúng. Trong các bước trình tự cơ bản để cho một nghiên cứu phát sinh sinh loài thì đánh giá sự phát sinh sinh loài cũng là một bước không thể bỏ qua. Sau đây là một số phương pháp được sử dụng hiện nay: Phương pháp Hà tiện tối đa (Maximum parsimony), một sự giả định cho rằng cây tiến hóa tốt nhất mổ tả tiến trình tiến hóa tốt nhất chính là cây mô tả được các loài ít thay đổi nhất tức là có ít đột biến nhất, cây vì thế có điểm thấp nhất (hà tiện) theo một tiêu chuẩn định sẵn. Phương pháp Khoảng cách (Distance method): Khác với phương pháp parsimony có mô hình tiến hóa là một hàm ẩn, thì phương pháp khoảng cách lại có mô hình tiến hóa là một hàm hiện. Trong phương pháp này từng cặp trình tự một sẽ được so sánh thẳng hàng cặp đôi và ứng với từng cặp, khoảng cách di truyền sẽ được tính toán. Do mô hình tiến hóa là một hàm hiện nên một trong số mô hình tiến hóa có thể được chọn để tính toán khoảng cách di truyền giữa từng cặp taxa từ đó cho ra một ma trận khoảng cách giữa tất cả các taxa. Và để có được cây tiến hóa, phương pháp phân rã hình ngôi sao thường được sử dụng ví dụ phương pháp neighbor-joining(liên kết cận kề). Do phương pháp neighbor-joining mà một trong những phương pháp nhanh nhất để dò tìm cây tiến hóa nên nó thường được sử dụng để phân tích khối dữ liệu lớn với nhiều taxa. Phương pháp Hợp lý cực đại (Maximum Likelihood) là phương pháp tiêu tốn nhiều thời gian nhất nhưng lại cho kết quả đáng tin cậy nhất. Mô hình tiến hóa dùng trong phương pháp này cũng là một hàm hiện. Ứng với mỗi mô hình tiến hóa được chọn, phương pháp này sẽ tính toán khả năng xác suất mà một cây tiến hóa có thể có từ chuỗi trình tự phân tích. Cây tiến hóa có xác suất cao nhất là cây cuối cùng được chọn. Chúng ta tập trung vào phương pháp ML, để hiểu được điều này chúng ta bắt đầu với những ví dụ cụ thể để mô hình hóa bài toán trên cây sinh loài. Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 32 4.3. Mô hình ước lượng hợp lý cực đại trên cây sinh loài Cho S1, S2, …, SN là một dãy mẫu DNA mà chúng ta có. Để đơn giản, giả thiết rằng mọi chuỗi trên có cùng chiều dài. Chúng ta muốn xác định những tham số của một cây sinh loài thông qua dãy mẫu trên và làm cực đại khả năng có thể xảy ra. Để giải bài toán này ta cần chỉ rõ một mô hình xác suất. Cho đơn giản, giả thiết “DNA” của chúng ta chỉ có hai trạng thái X và Y. Cạnh e được gán xác suất , có nghĩa là xác suất những thay thế (X ÙY) ngang qua e là (Hình 4). ep ep Hình 4: Mô tả xác suất chuyển đổi trạng thái của chuỗi “DNA” Phải chăng cạnh e được gán xác suất , có nghĩa xác suất của những mẫu liên quan thay thế ngang qua e, ví dụ XXYXYÙYXYXX được xác định rõ, và dễ dàng tính toán hàm Likelihood cho mẫu này: . ep 2 3(1 )e ep p− Qua bài toán trên có câu hỏi đặt ra như sau: Cái gì “hợp lý” mẫu trên? Có nghĩa là tìm kiếm mà nó làm cực đại xác suất của các mẫu trên. ep Mở rộng mô hình bài toán trên, mô hình mới của chúng ta sẽ gồm có một cây thông thường, nhưng ngoài ra các cạnh được gán những xác suất thay thế. Ví dụ ở đây, cây có 4 taxa. Những taxa này là sinh vật hoặc là gen, mỗi một taxa được mô tả bởi chuỗi DNA: Human: ATGGCTATTCTTATAGTACG Mouse: BATCGCTAGTCTTATATTACA Rate: TTCACTAGACCTGTGGTCCA Chicken: TTGACCAGACCTGTGGTCCG Hình 5: Cây sinh loài với các nút trong và xác suất chuyển đổi Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 33 Bây giờ chúng ta không biết trạng thái ở tại nút trong, đồng thời cũng không biết những tham số cạnh (Hình 5). 1 2 3 4 5 5 , , , , , ,e e e e e e ep p p p p p p 6 Hai hướng được đưa ra: 1. Cực đại qua những trạng thái của những nút bên trong. 2. Trung bình qua những trạng thái của những nút bên trong. Trong cả hai trường hợp, chúng ta đều làm cực đại những tham số qua cạnh. Trong hướng đầu tiên (trung bình, hoặc tổng những trạng thái những nút trong) chúng ta đang tìm kiếm “thích hợp nhất” đặt trên những cạnh của cây. Hướng này được gọi là cực đại khả năng cây sinh loài. Trong hướng này ML có lẽ là phương pháp suy diễn rộng rãi nhất được sử dụng hiện nay. Trong hướng thứ hai (làm cực đại qua những trạng thái của những nút trong) Chúng ta đang tìm kiếm “thích hợp nhất” những trạng thái tổ tiên. Hướng này được cực đại khả năng xảy ra ở tổ tiên (ancestral maximum likelihood -AML). Hướng thứ hai cũng phải sử dụng phương pháp ML bởi vì mục tiêu cuối cùng cũng phải là cực đại khả năng. 4.4. Mô hình tiến hóa Trong sinh vật học, quá trính tiến hóa là một quá trình phức tạp. Trong quá trình đó, các chuỗi gen phân kỳ từ cùng một tổ tiên. Nhưng vì sự đột biến và chia rẽ của sự đột biến đó làm tiến hóa cộng đồng bởi sự chọn lọc. Kết quả là sự thay đổi trạng thái của một nucleotide này thành một nucleotide khác ở những vị trí khác nhau. Trong việc tái cấu trúc cây sinh loài, chúng ta cần phải chấp nhận mô hình với một số giả định về quá trình cũng như trạng thái thay thế sau: - Mô hình đơn giản nhất là mô hình mà trong đó khả năng của bất kỳ nucleotide nào thay đổi thành bất kỳ nucleotide khác là bằng nhau. - Dự đoán khả năng rằng một nucleotide cụ thể ở một vị trí cụ thể sẽ thay đổi thành một nucleotide xác định khác trong một khoảng thời gian, cái chúng ta cần biết ở đây là tỷ lệ tức thời của sự thay đổi. Ma trận tỷ lệ (hoặc ma trận Q) là ma trận vuông ( ijQ q )= , với chỉ mục hàng và cột cho bởi { , , , }A C G T∑ = . Chúng ta cũng có thể sử dụng ký tự nhị phân hoặc 20 kí tự của amino axit cho tập . Ma trận tỷ lệ phải thỏa những yêu cầu sau: ∑ 0ijq ≥ với i j≠ 0ij j q ∈∑ =∑ cho tất cả i∈∑ , 0iiq < với tất cả i∈∑ Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 34 Ma trận tỷ lệ có được từ ý nghĩa từ tỷ lệ tức thời của đột biến. Từ ma trận tỷ lệ Q, chúng ta có thể tính được ma trận thay thế ( )tθ bởi hàm mũ theo công thức sau: 0 1( ) ! Qt i i i t e Q t i θ ∞ = = =∑ Phần tử của ( )tθ ở dòng i và cột j chính là xác suất mà sự thay đổi xảy ra trong một khoảng thời gian là t. i → →" j Mô hình đơn giản có một tham số và được biết là mô hình Jukes-Cantor với tỷ lệ chuyển đổi từ một nucleotide này đến một nucleotide khác là bằng nhau như sau: 3 3 3 3 A G C T A G Q C T α α α α α α α α α α α α α α α α −⎛ ⎞⎜ ⎟−⎜ ⎟= ⎜ ⎟−⎜ ⎟−⎝ ⎠ Và ma trận thay thế tương ứng là: 4 4 4 4 4 4 4 4 4 4 4 4 4 1 3 1 1 1 1 1 3 1 11( ) 4 1 1 1 3 1 1 1 1 1 3 t t t t t t t t t t t t t e e e e e e e e t e e e e e e e e α α α α α α α α α α α α α θ − − − − − − − − − − − − − − − − ⎛ ⎞+ − − −⎜ ⎟− + − −⎜ ⎟= ⎜ ⎟− − + −⎜ ⎟⎜ ⎟− − − +⎝ ⎠ 4 4 4 t t t α α α Ma trận ( )tθ thỏa: các phần tử của ma trận đều lớn hơn hoặc bằng 0 và nhỏ hơn hoặc bằng 1, tổng các phần tử trên một hàng bằng 1. Chúng ta cần xác định tính hợp lý mô hình trên. Giả sử chúng ta có G ở vị trí nào đó ở thời điểm t = 0, chúng ta hỏi rằng khả năng bao nhiêu ở đó vẫn là G vào thời điểm t (kí hiệu ), và tương tự như vậy khả năng là bao nhiêu nếu như A thay thế vào vị trí đó (kí hiệu ). Nếu tỉ lệ thay đổi là α trên đơn vị thời gian như mô hình Jukes - Cantor trên, thì: ( ) ( )GGP t ( ) ( )GAP t -4 ( ) 1 3( ) 4 4 tGGP t e α= + và -4( ) 1 1( ) 4 4 t GAP t e α= − Cũng theo mô hình Jukes-Cantor thì tất cả thay thế là như nhau, nên phát biểu chung là: -4 ( ) 1 3( ) 4 4 tiiP t e α= + và -4( ) 1 1( ) 4 4 t ijP t e α= − Ta thấy: Khi thì và , 0t → ( ) ( ) 1iiP t → ( ) ( ) 0ijP t → Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 35 Khi t thì →∞ ( ) 1( ) 4iiP t → và ( ) 1( ) 4ij P t → Điều này rõ ràng là phù hợp với thực tế của mô hình Jukes - Cantor, vì tại một thời điểm tức thời việc giử nguyên trạng thái với xác suất rõ ràng là bằng 1 và chuyển trạng thái này sang trạng thái khác là 0. Tương tự khi thời gian vô cùng lớn thì việc chuyển đổi trạng thái từ một nucleotide này sang một nucleotide khác hay giữ nguyên trạng thái là bằng nhau và bằng 0.25. Hiện nay, ngoài mô hình Jukes-Cantor còn có một số mô hình khác thường sử dụng như: Kimura-2, Kimura-3,…. Trong các mô hình này có sự khác nhau về tỉ lệ thay đổi trạng thái giữa các cơ sở. Khi sử dụng mô hình tiến hóa để tái cấu trúc cây, một là gán giá trị cụ thể cho tỉ lệ hoặc là ước lượng giá trị từ dữ liệu. Những mô hình này hoàn toàn giả định rằng các tốc độ là như nhau ở tất cả các vị trí. ML cố gắng suy ra một cây sinh loài bằng cách tìm ra cây mà cực đại khả năng đối với dữ liệu mẫu. Ví dụ: Dữ liệu mẫu ở đây là những chuỗi bằng nhau của nucleotides hoặc amino acids (chiều dài mỗi chuỗi N=32): TCAAAAATGGCTTTATTCGCTTAATGCCGTTA TCCGTGATGGATTTATTTCTGCAATGCCTGTC TTCGTGATGGATTTATTGTTGGTATGCCAGTC TTCGTGACGGGTTTATCTTGGCAATGCCGGTC Chúng ta bắt đầu với một mô hình tiến hóa cho bởi ma trận ( )tθ và một giả định một số hình dáng cây với chiều dài tương ứng. Có 15 khả năng cho các dạng cây có gốc với 4 taxa, một trong những cây đó là hình 6, trong đó các đỉnh ở lá tương ứng với các dữ liệu dóng theo cột được đánh dấu đánh đậm trên 4 chuỗi trên. Hình 6: Một trong những cây sinh loài 4 taxa Chúng ta không biết các nucleotide ở nút X và Y, nhưng có 4 khả năng xảy ra cho mỗi nút X và Y, vậy có có 16 trường hợp có thể xảy ra ở cây trên, một trong những trường hợp đó là hình 7. Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 36 Hình 7: Cây sinh loài với dữ liệu trên nút lá và các khả năng xảy ra ở các nút tổ tiên Xác xuất cho sự kiện mà mẫu quan sát A ở gốc là PA, bằng tần suất xuất hiện của A và thường được lấy bằng 0.25, giá trị này xác định được do thực nghiệm và độc lập với mô hình. Xác suất chuyển đổi từ A ở gốc đến G ở lá được tính toán từ ma trận ( )tθ và chiều dài của nhánh từ A đến G là ( )AGP . Vậy xác suất của cây là: 1 ( ) ( ) ( ) ( ) (. . . . . tree A AG AC AT TT TP P P P P P P= )T . Bởi vì có 16 trường hợp như vậy, xác suất của cây được tính bằng tổng các khả năng như sau: _ 1 2 tree i tree tree treeP P P P= + + +" 16 Đây chỉ là xác suất cho cây với dữ liệu quan sát ở một vị trí i được đánh dấu màu đậm ở các chuỗi trên. Khả năng của toàn bộ dữ liệu mẫu ở tất cả các vị trí bằng tích các khả năng cho mỗi một vị trí từ 1 đến N _ 1 N TREE tree i i P P = =∏ Áp dụng phương pháp ML là tìm xác suất các _tree iP hay rõ hơn là tìm xác suất chuyển đổi trạng thái trên các nhánh của cây dựa theo tham số của ma trận ( )tθ trên từng nhánh của cây sao cho xác suất đạt giá trị lớn nhất. TREEP Tuy nhiên việc giải bài toán trên là việc giải một hệ thống phương trình phi tuyến với nhiều ẩn số. Việc giải bằng tay là một điều không thể. Hướng giải quyết đưa ra là chọn một phương pháp toán học thích hợp kết hợp những ứng dụng của đại số máy tính hiện nay chúng ta có thể giải quyết bài toán trên với một số cây sinh loài nhỏ với một số mô hình chuyển đổi thường sử dụng.Với cách tiếp cận như vậy chúng ta có thể giải tìm nghiệm chính xác cho bài toán trên. Một trong những phương pháp đó là tìm thành phần bất biến trên cây sinh loài sẽ giới thiệu ở chương sau. Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 37 Chương 5. BẤT BIẾN TRÊN CÂY SINH LOÀI Theo chương trước chúng ta nhận thấy, giải bài toán cây sinh loài dẫn đến giải một bài toán cực đại một phương trình phi tuyến rất nhiều ẩn. Việc làm này khó khả thi ngay cả những cây sinh loài nhỏ. Người ta nhận thấy, đối với những mô hình thường sử dụng, trên cây sinh loài tồn tại những thành phần bất biến. Đối với một cây sinh loài cụ thể thành phần không đổi và không phụ thuộc vào mẫu dữ liệu quan sát. Từ những thành phần bất biến này, thay vì giải bài toán trên các tham số thì ta giải bài toán tương đương dựa trên các thành phần bất biến với sẽ đơn giản hơn. Trong chương này tập trung vào việc tìm tất cả thành phần bất biến. Cuối chương có một ví dụ về bất biến trên một cây sinh loài cụ thể. 5.1. Dẫn nhập Mô hình thống kê đại số mà chúng ta đang xét trên cây sinh loài là một ánh xạ có dạng: : d mf →^ ^ ở đây không gian tổng quát trên trường số phức, tuy nhiên các toạ độ thực tế 1,..., mf f của hàm số là các đa thức có hệ số hữu tỷ, có nghĩa là 1 1,..., [ ,..., ]m df f θ θ∈_ . Sử dụng phương pháp hợp lý cực đại, với các chuỗi dữ liệu quan sát ta có phương trình hợp lý tương ứng sau: mu m u ffL ...11= (1) trong đó là các số nguyên dương. 1,..., mu u Mục tiêu cuối cùng là làm cực đại hàm L trên bằng cách giải các phương trình đạo hàm riêng của nó: 0 0 m i i i i j u f f θ= ∂ =∂∑ với j = 1, .., d (2) Tuy nhiên, khi nghiên cứu bài toán trên cây sinh loài các đa thức 1,..., mf f là những phương trình phi tuyến với nhiều tham số 1, ..., dθ θ , cho nên việc giải bài toán trên với việc giải hệ phương trình (2) là việc làm khó khả thi. Một câu hỏi được đặt ra: các tập ảnh như thế nào khi các mff ,...,1 dθθ ,...,1 chạy trên miền xác định của nó? Nếu chúng ta xác định được tập các ảnh thì thay vì giải quyết bài toán (1) trên các tham số mff ,...,1 dθθ ,...,1 chúng ta chỉ xét bài toán (1) trên các tọa độ sẽ đơn giản hơn nhiều. mff ,...,1 Ví dụ: Chúng ta xét các ví dụ với ánh xạ 2 3:f →^ ^ Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 38 (i) Nếu thì tập ảnh là nghiệm phương trình ),,( 2121 2 1 θθθθθ=f 032 =− ff (ii) Nếu thì tập ảnh là nghiệm phương trình ),2,( 2221 2 1 θθθθ=f 04 31 2 2 =− fff (iii) Nếu thì tập ảnh là nghiệm phương trình ),,( 4221 4 2 5 121 5 1 θθθθθθθθ +++=f 0)()()(2 20231 5 132 4 321 11 =−+−−+−+ fffffffff Qua các ví dụ trên có nhận xét sau: Các đa thức ở dạng những đơn thức thì tập ảnh là nghiệm của những phương trình đơn giản, ngược lại thì là những phương trình bậc rất lớn. mff ,...,1 Quay trở lại bài toán trên cây sinh loài chúng ta quan tâm, những đa thức đối với từng cây cụ thể có tính chất riêng của chúng. Tập hợp ảnh của những tức là tìm tất cả các phương trình quan hệ của chúng. Tập các phương trình như thế được gọi là thành phần bất biến của cây sinh loài. mff ,...,1 mff ,...,1 Mục tiêu của phần này giới thiệu một phép biến đổi có tên là phép biến đổi Fourier để tìm ra tất cả các bất biến trên một cây sinh loài cụ thể. 5.2. Mô hình xác suất trên cây sinh loài 5.2.1. Mô hình bài toán cây sinh loài Cho T là cây có gốc với n lá. Đặt V(T) là tập các nút của T. Với mỗi một , chúng ta kí hiệu biến Xv, mỗi biến này mang 1 trong k giá trị. Trong sinh học, k hầu như có các giá trị 2, 4 và 20. Kí hiệu P(Xv = i) cho xác suất Xv mang trạng thái i. ( )v V T∈ Mối quan hệ giữa các biến ngẫu nhiên Xv được xác định bởi cấu trúc của cây. Đặt π là phân bố của biến Xr tại nút gốc r. Với mỗi một nút ( ) \ { }v V T r∈ , đặt a(v) là nút cha duy nhất của v. Sự chuyển trạng thái từ a(v) đến v được cho bởi ma trận xác suất chuyển đổi A(v) có kích cỡ k k× . Và xác suất phân bố ở mỗi một nút được tính toán đệ quy như sau: ∑ = === k i va v ijv iXPAjXP 1 )( )( )(.)( Công thức này được suy ra từ phân bố trên tất cả biến ngẫu nhiên Xv. Chúng ta gán nhãn các nút lá cho T bởi 1, 2, …, n và ta có xác suất phân bố các biến tại các lá: ),...,,( 2211..21 nniii iXiXiXPp n ==== Trong các ứng dụng sinh học, người ta ước lượng có kn khả năng từ n chuỗi bằng nhau trên k kí tự. Mục đích chúng ta là dựa vào n chuỗi bằng nhau đó, xác định hình dáng cây sinh loài ở quá khứ mà khả năng xảy ra lớn nhất, nói cách khác là tái cấu trúc cây sinh loài. Vậy đầu vào bài toán chúng ta chỉ có n mẫu dữ liệu, tức là n chuỗi DNA tương ứng, các phân bố gốc π và ma trận chuyển trạng thái A(v) là chưa biết. Tuy nhiên để đơn giản cho các bài toán, người ta đưa ra các mô hình đơn giản Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 39 gần với thực tế thường sử dụng là: Phân phối π là phân phối đều và ma trận chuyển trạng thái A(v) được sử dụng là mô hình Jukes - Cantor hay Kimura 2 và 3 trạng thái. Với các giả thiết trên, các bất biến của cây sinh loài của mô hình là một đa thức dựa trên các khả năng ở lá là và triệt tiêu với mọi sự chọn lựa tham số của mô hình. Tập các đa thức là iđêan nguyên tố trên vành đa thức với các biến chưa biết . Mục tiêu chúng ta là tìm các iđêan này. niiip ..21 niiip ...21 5.2.2. Nhóm Abel và sự liên hệ với các ma trận chuyển đổi Ở mô hình Neyman trên 2 kí tự (k = 2) là mô hình với ma trận chuyển đổi ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − −= vv vvv aa aa A 1 1)( với av là xác suất được tạo ra bởi sự chuyển đổi giữa các trạng thái dọc theo cạnh từ a(v) đến v. Mô hình Kimura với 3 tham số với k = 4 kí tự (đối với chuỗi DNA) có ma trận chuyển đổi sau: ⎟⎟ ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜⎜ ⎜ ⎝ ⎛ −−− −−− −−− −−− = vvvvvv vvvvvv vvvvvv vvvvvv v cbaabc acbacb bccbaa cbacba A 1 1 1 1 )( Mô hình Kimura 2 tham số được định nghĩa như ma trận trên với v vb c= . Tương tự, mô hình Jukes- Cantor với 4 kí tự, ma trận trên với . v va b c= = v Chìa khóa đối với trạng thái các biến ngẫu nhiên vX là nhóm hữu hạn Abel (ví dụ hoặc với phép toán cộng trên các tọa độ và mod cho 2). Giả sử rằng, chúng ta xem các cơ sở {A, G, C, T} như là các phần tử của nhóm Abel, với phép toán được định nghĩa với bảng cộng sau: 2 {0,1}=Z 2 2 {(0,0),(0,1),(1,0),(1,1)}⊕ =Z Z AG C T GA T C C TA G T CG A T C G A TCGA+ Nhóm trên đẳng cấu với nhóm 2 2⊕Z Z , tương ứng với AÙ(0,0), GÙ(0,1), CÙ(1,0), TÙ(1,1). Từ đó dễ thấy rằng, ma trận ở mô hình Kimura 3 tham số có tính chất tương ứng với từng cặp cơ sở . Mặt khác, chúng ta thấy rằng sự chuyển trạng thái từ đến ( , )i jg g ig jg chỉ phụ thuộc vào hiệu ig g j− . Điều đó cũng đúng khi chúng ta xem xét Phương pháp đại số cho bài toán ước lượng hợp lý cực đại - Áp dụng trên cây sinh loài nhỏ Bùi Văn Đồng Trang 40 đối với mô hình Jukes-Cantor. Vì thế, ma trận chuyển đổi của những mô hình đang quan tâm không gì khác hơn những ma trận trên nhóm chúng ta đang xét. 5.3. Biến đổi Fourier Đặt là nhóm hữu hạn Abel với phép toán được viết như là +. Đặt G { :z z= ∈ =Y C 1} như là vòng tròn đơn vị trên mặt phẳng số phức. Chúng ta thấy rằngYlà một nhóm Abel với phép toán nhân thông thường của số phức. Những đặc trưng của là những đồng cấu nhóm từ vào Y. Nghĩa là G G :χ →G Y là một đặc trưng nếu 1 2 1 2