Khóa luận Giải pháp sắp xếp chuỗi gen sử dụng cây hậu tố

J. Watson và F.Crick làm việc tại trường đại học Cambridge đã xây dựng thành công mô hình không gian của phân tử DNA(deoxyribonucleic acid), đánh dấu một bước ngoặt quan trọng trong sự phát triển của sinh học phân tử, theo mô hình này DNA là một đại phân tử sinh học có cấu trúc như một chuỗi xoắn kép gồm hai mạch đơn, mỗi mạch đơn là một chuỗi nucleotide. Mỗi nucleotide gồm nhóm phosphate, đường desoxyribose và một trong bốn thành phần lần lượt được biểu thị bởi các chữ cái A, C, G và T. Hai mạch đơn kết hợp với nhau nhờ các liên kết hydro hình thành giữa các thành phần bổ sung nằm trên hai mạch. A bổ sung cho T, C bổ sung cho G. Bảng 2: các nucleotide[23] DNA Adenine Guanine Cytosine Thymine A G C T/U RNA Adenine Guanine Cytosine Uracil Mỗi mạch đơn là một trình tự có định hướng với một đầu là đầu 5’ phosphate tự do còn đầu kia là đầu 3’ hydroxyl tự do (hướng quy ước là 5’-> 3’) . Hướng của hai mạch đơn trong chuỗi xoắn kép ngược nhau, người ta gọi chúng là hai mạch đối song song. Từ định nghĩa trên nảy sinh ra hai khái niệm cơ bản: Mỗi mạch đơn là một trình tự những thành phần khác nhau nên mỗi mạch đơn mang thông tin khác với mạch còn lại. Hai mạch đơn liên kết với nhau bởi một quan hệ bổ sung, chính quan hệ này giải thích được cấu trúc chặt chẽ của phân tử DNA và đặc biệt là cách tự sao chép để tạo ra hai phân tử con từ một phân tử mẹ. Hình 3: Cấu trúc xoắn kép DNA DNA đóng vai trò cơ bản trong quy trình của sự sống dưới ở 2 phương diện. Trước hết nó chứa khuôn mẫu cho sự tổng hợp proteins, điều này thiết yếu với bất kì sinh vật sống nào. Mặc dù có khá nhiều loại protein khá rộng lớn nhưng thành phần chung tạo nên protein chính là các axit amin. Mỗi axit amin trong 20 axit amin đó được mã hóa bởi một hoặc nhiều bộ ba nucleotides tạo nên DNA, cứ 3 nucleotide sẽ mã hóa 1 loại amino axid Dựa trên bảng dịch mã, một xâu DNA tuyến tính được dịch sang một xâu axit amin tuyến tính. Dưới đây là một ví dụ: ... Hình 4: Sơ đồ mã hóa amino axid[23] Vai trò thiết yếu đối với cuộc sống thứ 2 của DNA là bảo quản và truyền đạt thông tin về các loại protein trong cơ thể, nói cách khác nó lưu trữ thông tin di truyền của sinh vật. Trong cấu trúc của DNA, mỗi cặp sợi bổ sung cụ thể là (A-T và G-C) tạo thành một xoắn kép. Vì thế mỗi sợi mang toàn bộ thông tin và Bộ máy hóa sinh đảm bảo rằng thông tin có thể sao lại hết lần này đến lần khác ngay cả khi phân tử gốc đã biến mất từ lâu. Trong quá trình sao chép, sự thay đổi hay còn gọi là đột biến có thể xảy ra đối với chuỗi DNA. Phân loại đột biến rất quan trọng trong việc so sánh chuỗi dựa trên sự thay đổi, sự chèn nucleotide vào trong chuỗi và sự loại bỏ nucleotide ra khỏi chuỗi. Hoạt động sơ cấp được cho phép trong định nghĩa chuỗi tương tự là lựa chọn để các chuỗi trở nên tương ứng với nhau. Để hình dung mối quan hệ giữa 2 chuỗi tương tự ta có thể tham khảo ví dụ của sự sắp xếp sau: V-LSPADKTNVKAAWGKVGAHAGEYGAEALERMFLSFPTTKTYFPHF-DL HAHU VHLTPEEKSAVTALWGKV--NVDEVGGEALGRLLVVYPWTQRFFESFGDL HBHU SH-----GSAQVKGHGKKVADALTNAVAHVDDMPNALSALSDLHAHKLRV HAHU STPDAVMGNPKVKAHGKKVLGAFSDGLAHLDNLKGTFATLSELHCDKLHV HBHU DPVNFKLLSHCLLVTLAAHLPAEFTPAVHASLDKFLASVSTVLTSKYR HAHU DPENFRLLGNVLVCVLAHHFGKEFTPPVQAAYQKVVAGVANALAHKYH HBHU Hai chuỗi axit amin được so sánh ở đây là chuỗi alpha và beta của hemoglobin người (lần lượt viết tắt là HAHUvà HBHU). Với chuỗi dài xấp xỉ 150 axit amin, mỗi khối dòng chứa một phần của chuỗi thứ nhất ở dòng trên và của chuỗi thứ 2 ở dòng dưới. Phần còn lại ở trên mỗi chuỗi trong một khối là tương đương. Một số được giữ lại (axit amin trong cột giống hệt nhau), và một số trao đổi và một phần của chuỗi bị xóa khỏi một chuỗi tương đương với việc chèn thêm vào chuỗi kia. Chèn và xóa được chỉ định bởi một cặp chữ cái với dấu gạch ngang. Sự sắp xếp là một hoạt động cần thiết để biến đổi một chuỗi thành chuỗi khác sử dụng cùng một hành động tương tự như sự tiến hóa. DNA có khả năng tự nhân đôi, quá trình này xảy ra chủ yếu bên trong nhân tế bào, tại các nhiễm sắc thể ở kì trung gian giữa hai lần phân bào. Trong giai đoạn này, chuỗi DNA tự tách ra thành hai chuỗi đơn, các chuỗi đơn kết hợp với nucleotide trong môi trường nội bào để tạo thành một chuỗi đơn mới theo nguyên tắc bổ sung. Sau đó, các chuỗi đơn ban đầu và chuỗi mới xoắn lại với nhau lại với nhau tạo thành 2 phân tử DNA giống hệt phân tử DNA gốc ban đầu. Trong quá trình DNA tách ra làm đôi, nó có thể bị tác động bởi các tác nhân hóa, lý làm thay đổi (mất đi, chèn thêm nucleotide) trên các chuỗi đơn. Như vậy DNA mới sinh sẽ không còn giống như DNA gốc, và những thay đổi này sẽ được biểu hiện ngay thành các tính trạng sinh vật từ thế hệ này sang thế hệ khác. Điều này giải thích cho quá trình tiến hóa không ngừng của sinh vật. Có một giả thuyết rất nổi tiếng đã được đưa ra là:”Tất cả các sinh vật đều có cùng một tổ tiên!”, tức là đều có chung nguồn gốc di truyền. Cách đây vài thập kỷ, người ta chỉ có thể hình dung ra DNA trong tưởng tượng tuy nhiên giờ đây ta có thể quan sát DNA một cách rất trực quan như hình 5: Hình 5: Hình ảnh về DNA[3] RNA Giống như DNA, RNA cũng có cấu trúc đa phân mà đơn phân là 4 loại nucleotide, tuy nhiên trong RNA nucleotide loại T (pyrimidine thymine) được thay thế bằng U (uracil). RNA tồn tại ở dạng chuỗi đơn và được phân chia làm 3 loại chính dựa trên chức năng của chúng: mRNA (RNA thông tin): là một mạch sao chép nguyên từ một mạch đơn của DNA trong đó T được thay bằng U và làm nhiệm vụ truyền đạt thông tin cấu trúc protein được tổng hợp. rRNA (RNA riboxom): là thành phần cấu tạo nên riboxom. tRNA (RNA vận chuyển): có chức năng vận chuyển amino axid tương ứng đến nơi tổng hợp protein. snRNA: có chức năng hỗ trợ việc ghép mã mRNA. gRNA: sử dụng để điều khiển việc thay đổi mRNA . RNA có thể liên kết với một dải đơn của một phân tử DNA, bằng cách thay T bằng U, và các phân tử kiểu này có vai trò quan trọng trong các quá trình sống và công nghệ sinh học.[1] C-G-A-T-T-G-C-A-A-C-G-A-T-G-C DNA | | | | | | | | | | | | | | | G-C-U-A-A-C-G-U-U-G-C-U-A-C-G RNA Hình 6: Hình ảnh về RNA[10] Protein Protein là một đại phân tử sinh học được hình thành từ 1 hay nhiều chuỗi polypeptide sắp xếp theo một thứ tự đặc biệt, thứ tự này được xác định bởi dãy cơ sở (peptide là một chuỗi nối tiếp nhiều axit amin với số lượng ít hơn 30, với số lượng axit amin lớn hơn chuỗi được gọi là polypeptide) được hình thành từ 20 loại axit amin khác nhau lần lượt được biểu thị bằng 20 kí tự khác nhau trong bảng chữ cái. Từ “protein” dùng để chỉ một cấu trúc phức tạp trong không gian chứ không đơn thuần chỉ là một trình tự axit amin. Các nucleotides trong gene mã hóa cho protein. Các proteins cần thiết cho cấu trúc, chức năng và điều chỉnh tế bào, mô và tổ chức, mỗi protein có một vai trò đặc biệt. Bảng 3: 20 axit amin[23] STT One-letter code Three-letter-code Name 1 A Ala Alanine 2 C Cys Cysteine 3 D Asp Aspartic Acid 4 E Glu Glutamic Acid 5 F Phe Phenylalanine 6 G Gly Glycine 7 H His Histidine 8 I Ile Isoleucine 9 K Lys Lysine 10 L Leu Leucine 11 M Met Methionine 12 N Asn Asparagine 13 P Pro Proline 14 Q Gln Glutamine 15 R Arg Arginine 16 S Ser Serine 17 T Thr Threonine 18 V Val Valine 19 W Trp Tryptophan 20 Y Tyr Tyrosine Bảng 4: thứ tự các nucleotide trong axit amin[3] First Position Second Position Third Position T C A G T Phe Ser Tyr Cys T Phe Ser Tyr Cys C Leu Ser Stop Stop A Leu Ser Stop Trp G C Leu Pro His Arg T Leu Pro His Arg C Leu Pro Gln Arg A Leu Pro Gln Arg G A Ile Thr Asn Ser T Ile Thr Asn Ser C Ile Thr Lys Arg A Met Thr Lys Arg G G Val Ala Asp Gly T Val Ala Asp Gly C Val Ala Glu Gly A Val Ala Glu Gly G Mỗi amino axid bao gồm một nguyên tử Cacbon trung tâm (Cα), một nhóm amino (NH2), một nguyên tử Hidro, một nhóm COOH và một gốc Ri liên kết trực tiếp với nguyên tử Cα. Gốc Ri là đặc trưng cho từng loại amino axid khác nhau. Hình 7: Liên kết peptit giữa các amino axid[13] Protein bao gồm bốn mức độ tổ chức: Cấu trúc bậc 1 là trình tự sắp xếp các axit amin trong chuỗi polypeptid, cấu trúc bậc 2 phát sinh từ sự uốn các thành phần của chuỗi polypeptid thành những cấu trúc đều đặn trong không gian( dạng xoắn α(alpha helix) hay lớp mỏng β(Beta sheets)). Cấu trúc bậc 3 quy định sự kết hợp các chuỗi xoắn hay lớp mỏng đó thành hình dạng ba chiều trong không gian. Cấu trúc bậc 4 là sự tổ chức nhiều chuỗi polypeptid thành một phân tử protein. DNA, RNA và quá trình tổng hợp protein Tổng hợp protein là quá trình tạo ra proteins dựa trên thông tin được mã hóa trong gen( là các đoạn mã đặc biệt của DNA có chức năng điều khiển cấu trúc và hoạt động của tế bào, là đơn vị chức năng của sự di truyền) gồm ba giai đoạn chính: (1) Transcription (phiên mã) (2) Splicing (ghép mã) (3) Translation (dịch mã)[1] có thể được mô tả như hình dưới Hình 8: Quá trình tổng hợp protein[1] Chương 2: Sắp xếp hai chuỗi sinh học Các chuỗi sinh học có thể được hiểu là dạng biểu diễn tuyến tính của các đại phân tử sinh học như axid nucleic (DNA, RNA) hay protein…Trong một chuỗi sinh học, các đơn phân của đại phân tử sinh học được ký hiệu bởi các chữ cái đại diện tương ứng. Như vậy, một chuỗi sinh học sẽ bao gồm một tập các chữ cái đại diện cho các phần tử đơn phân tương ứng. Chẳng hạn, đối với chuỗi sinh học của axid nucleic sẽ là một dãy các kí tự A,T,G,C xen kẽ lẫn nhau, đối với chuỗi sinh học của protein sẽ bao gồm 20 kí tự tương ứng với 20 loại axid amin khác nhau. Do DNA hay Protein đều có cấu trúc đa phân được cấu thành từ một số lượng rất lớn đơn phân, chính vì vậy các chuỗi sinh học cũng có kích thước rất lớn, thông thường cỡ hàng triệu phần tử. Đây chính là nguồn gốc của sự bùng nổ dữ liệu trong Tin Sinh Học. Hầu hết mọi người tin rằng sở dĩ các sinh vật có cùng tổ tiên trở nên khác nhau là do trải qua 1 quá trình tiến hóa. Tiến hóa là một quá trình mà kết quả của nó là sự thay đổi các đặc tính được thừa kế qua nhiều thế hệ. Bởi các thay đổi có khả năng được thừa kế nên các nhà khoa học tổng kết rằng phải có một vài thay đổi trong gen của các sinh vật này tương ứng với mỗi sự thay đổi tiến hóa. Trong thực nghiệm các nhà sinh học đã so sánh gen của hai sinh vật để hiểu rõ hơn về mối quan hệ di truyền và liên kết tiến hóa giữa 2 sinh vật. Điều này cho thấy rằng hai loài càng gần nhau càng có nhiều cặp gen chung. Ví dụ dưới đây về loài người và chuột sẽ cho ta một vài hình dung về điều này Bảng 5: Số cặp gen chung giữa người và chuột[22] Bởi lý do trên chúng ta phải tìm ra phương pháp so sánh bộ gen của hai loài và rút ra tất cả các cặp gen chung được bảo tồn của chúng. Để giải quyết vấn đề này chúng ta sẽ làm quen với việc so sánh 2 chuỗi sinh học. Giới thiệu về so sánh 2 chuỗi sinh học Việc nghiên cứu dựa trên phân tích các chuỗi sinh học đơn lẻ có nhiều hạn chế trong việc đưa ra các kết quả theo ý muốn đặc biệt không thể đánh giá mối tương quan giữa các chuỗi sinh học hay chính là mối quan hệ di truyền giữa loài này với loài khác. Chính vì vậy người ta đã đi đến phân tích trên hai hay nhiều chuối sinh học khác nhau. Có nhiều lý giải cho việc so sánh chuỗi này. Trước hết, thuyết tiến hoá cho rằng các chuỗi gen có thể đều xuất phát từ các chuỗi di truyền ban đầu. Do đó việc theo dõi lịch sử tiến hoá của những đột biến và các biến đổi tiến hoá khác là một công việc hết sức thú vị. Sự so sánh về các chuỗi sinh học trong trường hợp này được hiểu là sự so sánh dựa vào tiêu chí tiến hoá. Ví dụ, số lượng đột biến, chèn và xoá các thành phần cần thiết để biến đổi một chuỗi DNA thành 1 chuỗi khác là sự tính toán phản ánh mối quan hệ tiến hoá.  Nói cách khác, 1 phép so sánh có thể thực tế hơn khi mà người ta không nhằm vào 1 sự tái lập cụ thể quá trình tiến hoá của các sự kiện mà tập trung vào các đoạn xác định của nguồn gốc chung có thể lần lượt trùng hợp với các đoạn cấu trúc hoặc chức năng tương tự. Theo quan điểm này, các đặc tính vật lý của axit amin đóng vai trò quan trọng hơn trong nghiên cứu về tiến hoá. Biểu đồ điểm Biểu đồ điểm là phương pháp lâu đời nhất để so sánh các chuỗi (do Maizel và Lenk đưa ra). Một biểu đồ điểm là sự biểu thị các điểm tương đồng giữa 2 chuỗi sinh học. Biểu đồ gồm 2 trục, mỗi trục biểu thị một trong 2 chuỗi được so sánh. Một cửa sổ có độ dài cố định, theo một tiêu chí, là khi cho 2 chuỗi trượt trên các trục tương ứng, tại mỗi vị trí của chúng, các điểm(đường chéo ngắn) được tạo ra khi có sự giống nhau giữa các cặp phần tử trong 2 chuỗi tại vị trí đó. Vì vậy, nếu hai chuỗi giống nhau trên toàn bộ chiều dài, thì một đường chéo chính sẽ được tạo ra gồm tập hợp các điểm. Ngược lại, nếu 2 chuỗi chỉ tương tự nhau ở các đoạn bộ phận, thì một đường chéo kéo dài sẽ được tạo ra trên biểu đồ. Hình 10 đưa ra một ví dụ về biểu đồ điểm. Ở đây, chuỗi alpha Hemoglobin được so sánh với chuỗi beta hemoglobin ở người. Ở đây độ dài của ô được đặt là 31, các điểm tương xứng và không tương xứng được gán các giá trị tương ứng +5 đến -4. Các giá trị màu xám của biểu đồ sắp xếp theo sự tương đồng của 2 ô. Chúng ta có thể dễ dàng nhận thấy một đường chéo chính ở giữa được thể hiện rõ nét và các đường chéo mờ ở 2 bên, điều này cho thấy 2 chuỗi được so sánh là tương tự nhau trên toàn bộ chiều dài. Hình 10: Biểu đồ điểm của 2 chuỗi DNA: chuỗi alpha của hemoglobin ở người được quy định là trục hoành, còn trục tung là chuỗi beta . Hình được vẽ bởi công cụ dotter.[23]  Biểu đồ điểm là một phương pháp hữu hiệu để so sánh 2 chuỗi. Tuy nhiên đó không phải là phương pháp phân tích hiệu quả nhất. Dựa vào biểu đồ điểm chúng ta có thể đánh giá ngay 2 chuỗi được so sánh là tương tự trên toàn bộ chiều dài hay chỉ tương tự trên một đoạn đoạn cục bộ. Sự tương tự cục bộ có thể hiểu là sự tồn tại của các vùng trên hai chuỗi có sự tương tự nhau. Biểu đồ điểm sẽ chỉ ra tất cả các tương tự cục bổ tương ứng với những đường chéo nhỏ song song nhau ở hai bên đường chéo chính. Phương pháp này thể hiện một cách rất trực quan và dễ hiểu, thích hợp với những đánh giá mang tính chất định tính. Sắp xếp 2 chuỗi sinh học Sắp xếp chuỗi là sự sắp xếp một chuỗi ở trên một chuỗi khác nơi các phần tử nằm tại một vị trí được coi là có cùng nguồn gốc tiến hoá. Nếu các phần tử giống nhau cùng xuất hiện ở cả 2 chuỗi thì vị trí này được bảo tồn trong tiến hoá. Nếu các phần tử đó khác nhau, người ta cho rằng 2 chuỗi đó xuất phát từ cùng 1 phần tử ban đầu. Các chuỗi tương đồng có thể khác nhau về độ dài được giải thích bằng các thao tác chèn hay xoá trong các chuỗi sinh học. Thay vì đánh giá 2 chuỗi tương tự nhau một cách định tính như trong biểu đồ điểm chỉ ra, chúng ta có đánh giá một cách định lượng sự tương tự đó. Nói cách khác là chúng ta sẽ lượng tử hóa sự giống nhau giữa hai chuỗi nhằm đưa ra một kết quả đánh giá chính xác nhất. Để làm được việc này biểu đồ điểm rõ ràng không thể là một giải pháp thích hợp, thay vào đó chúng ta phải tính được giá trị biểu diễn độ giống nhau giữa 2 chuỗi đó dựa trên giá trị điểm số của từng cặp phần tử tương ứng nhau trong sắp xếp. Đối với sắp xếp các chuỗi DNA, chúng ta phải gán một giá trị điểm số cho tất cả các cặp có thể {A,A}, {T,T}, {A,G}, {C,G}… Còn trong trường hợp các chuỗi axid amin, có tới 20 loại axid amin khác nhau, chính vì vậy người ta thường biểu diễn điểm số của tất cả các cặp dưới dạng một ma trận điểm số. Đây là một ma trận 2 chiều trong đó mỗi một phần tử biểu diễn điểm số của một cặp 2 phần tử axid admin tương ứng. Giá trị điểm số có thể nhận các giá trị được qui định dựa trên những kinh nghiệm trong sinh học. Một điều dĩ nhiên là các cặp giống nhau {A,A},{M,M} sẽ có điểm số cao hơn so với các cặp không giống nhau {A,M},{K,V}…Có một số ma trận loại này đã được lập sẵn ra để sử dụng chung, có thể kể đến như: DayHoff, hay bộ ma trận trọng số BLOSUM… Dựa vào đây ta có thể tính điểm số cho sự giống nhau trên toàn chuỗi sinh học, đó là nền tảng cho các phương pháp so sánh chuỗi về sau này. Một trong những phương pháp đó gọi là sắp xếp chuỗi sinh học. Sắp xếp chuỗi sử dụng các phương pháp đối sánh 2 chuối sinh học với nhau để thể hiện sự tương đồng giữa chúng. Bằng cách gán thêm các khoảng trống(được biểu diễn bởi dấu ‘-‘) vào các vị trí nhất định trên mỗi chuỗi để tạo ra hai chuỗi được sắp xếp có sự tương ứng các cặp phần tử khác với ban đầu mà vẫn giữ được thứ tự các phần tử trên từng chuỗi. Chẳng hạn với hai chuỗi ban đầu BANANA và ANANAS ta có thể đưa ra một sắp xếp như sau: BANANA- -ANANAS Ký hiệu “-” dùng để biểu thị thao tác chèn hoặc xoá bởi vì chèn trong chuỗi này tương đương một sự bỏ trống trong một chuỗi khác, hiện tượng này gọi là “indel”, viết tắt của ‘insert’ và ‘delete’. Ta xét một ví dụ khác về sắp xếp chuỗi như dưới đây: B--KETBALL BASKE-BAL- Trong sắp xếp trên chúng ta cũng chỉ có thể đánh giá một cách định tính về sự tương đồng giữa 2 chuỗi, để có thể đưa ra đánh giá chính xác người ta đã sử dụng. phương pháp gán giá trị (điểm số) cho các cặp chữ được đối chiếu. Chẳng hạn đối với cặp chữ giống nhau sẽ gán cho một điểm số là một số dương (tùy theo từng cặp chữ cái), với các cặp không khớp nhau sẽ được gán điểm số âm, còn lại các cặp có một vế bị thiếu (indel) sẽ được cho một giá trị điểm số thấp. Sau đó tính tổng điểm số trên các cặp chữ cái sẽ cho ta điểm số của phương pháp trên. Việc tìm ra cách sắp để có thể đạt được điểm số tối đa sẽ là phương án sắp xếp tối ưu nhất. Như vậy, có thể thấy rằng chúng ta đã đưa bài toán trong sinh học về một bài toán có thể biểu diễn trong tin học. Bằng việc xây dựng các hàm tính điểm số, hàm trừ điểm số cho các vị trí trống và sử dụng các thuật toán trong tin học ta có thể giải quyết được bài toán trên. Một phương pháp cơ bản, dễ thấy là quét hết các khả năng để tìm trường hợp cho kết quả tốt nhất. Cách này chỉ dùng được khi kích thước các chuỗi là nhỏ. Phương pháp qui hoạch động sẽ cho kết quả khả quan hơn về thời gian tính toán và đã được áp dụng trong hầu hết trong các bài toán sắp xếp chuỗi sinh học. Tuy nhiên, một vấn đề rất đáng cân nhắc là việc xây dựng các hàm khoảng cách, hàm khoảng trống như thế nào là phù hợp nhất. Hàm khoảng cách sẽ xác định các giá trị âm của các chữ cặp khác nhau hoặc chỗ trống, sau đó cực tiểu hóa khoảng cách này. Hàm tương tự sẽ gán giá trị cao cho các cặp giống nhau và giá trị thấp của cho chỗ trống, sau đó cực tiểu hóa điểm. Việc này hoàn toàn là do cân nhắc dựa trên sự hiểu biết tổng quát về thế giới sinh học phân tử cũng như trên cơ sở toán học. Năm 1981, Smith và Waterman lần đầu tiên đã áp dụng thuật toán trên toàn chiều dài của chuỗi sinh học và đã thu được kết quả rất khả quan. Về sau này, phương pháp trên không ngừng được cải tiến về nội dung và cách thức thực hiện. Về vấn đề cân nhắc cho điểm số đối với các cặp chữ khớp nhau hay các cặp chữ không khớp nhau, chúng ta có thể dễ dàng thấy được sự khác biệt ngay trong bản thân chúng. Chẳng hạn, nếu một cặp Tryptophanes (T) giống nhau sẽ được đánh giá là quan trọng hơn so với cặp Analines (A). Đối với các axid amin, một ma trận điểm đối xứng sẽ được tạo ra để gán giá trị điểm số cho từng cặp tương ứng. Một trong những ma trận rất nổi tiếng là ma trận DayHoff ra đời năm 1978 mang tên người xây dựng nên nó. Một số ma trận khác được xây dựng gần đây hơn như các ma trận BLOSUM hay ma trận của Gonnet. Vấn đề còn lại lúc này là phải tìm ra một hàm đánh điểm số phù hợp cho các khoảng trống trên chuỗi. Trong thuật toán của Smith và Waterman đã không có một yêu cầu nào cụ thể trên các dãy khoảng trống có chiều dài khác nhau, nghĩa là nó coi một dãy các khoảng trống liên tiếp như là tổng của các khoảng trống đơn lẻ. Điều này có thể gây ra sự không chính xác trong cấu trúc tổng thể bởi người ta quan niệm rằng sự xuất hiện khoảng trống đầu tiên sẽ khác với các khoảng trống liền sau nó. Quan điểm này phù hợp với các qui luật tự nhiên và sinh học và đã được các nhà nghiên cứu xây dựng thành công thức dạng g(k) = a+b(k). Trong đó g là hàm khoảng trống, k là số khoảng trống kế tiếp, a và b là các hằng số. Việc xây dựng được các hàm khoảng trống và ma trận cặp điểm cùng với phương pháp qui hoạch động đã làm nền tảng để giải quyết bài toán sắp xếp chuỗi. Có rất nhiều thuật toán cài đặt cụ thể, ta sẽ làm quen với thuật toán sắp xếp loại I. Thuật toán sắp xếp 2 chuỗi sinh học Giới thiệu về thuật toán sắp xếp loại I Chúng ta sẽ làm quen với thuật toán dùng để so sánh 2 chuỗi, dựa vào quy hoạch động. Ngoài ra, còn đề cập tới các hàm trống, hàm tính điểm và khoảng cách. Đó là cách tiếp cận đối với sắp xếp 2 chuỗi sinh học bằng phương pháp sắp xếp loại 1, phương pháp này áp dụng kỹ thuật qui hoạch động một khá cách hiệu quả. Xét ví dụ sắp xếp 2 chuỗi sau: "RDISLVKNAGI" và "RNILVSDAKNVGI" R D I S L V - - - K N A G I  R N I - L V S D A K N V G I Các dấu “-” biểu thị khoảng trống trong chuỗi (phần chèn hay xoá ). Trật tự xuất hiện chữ cái trong sắp xếp giống nhau ở các chuỗi tương tự. Hai chữ cái in trong cùng cột gọi là cặp tương ứng. Một loại sắp xếp khác có thể thay thế cho 3 cặp đầu có thể là . Sắp xếp loại I không cho phép các chỗ trống liền kề trong chuỗi đối diện. Sắp xếp loại I có thể được biểu thị như 1 chuỗi các cặp chữ cái liên tục không ảnh hưởng đến trật tự của các chữ cái trong chuỗi. Ít nhất một trong số các chuỗi, không có chữ cái nào bị bỏ sót khi chuyển sang cặp tiếp theo. Ví dụ trên đây được coi là sắp xếp loại I được biểu thị bằng việc liệt kê các cặp số liệu đối với các cặp chữ cái. ( (1,1), (2,2), (3,3), (5,4), (6,5), (7,9), (8,10), (9,11), (10,12), (11,13)) . Trong ví dụ (5,4) có nghĩa là chữ cái thứ 5 của chuỗi đầu tiên tương ứng với chữ cái thứ 4 của chuỗi thứ 2. Cách trình bày này sẽ được chọn để xác định sắp xếp I. Trong 20 axit amin, các cặp chữ cái giống nhau không xảy ra thường xuyên như ở DNA. Vì vậy, sắp xếp trọng số đã được phát triển dựa trên thuộc tính giá trị của mỗi cặp axit amin phù hợp với nhau. Sắp xếp này do M. Dayhoff đưa ra dựa trên tần suất trao đổi giữa các axit amin. Ma trận thuộc tính 20 x 20 này xuất phát từ các giá trị dương khác nhau (sắp xếp từ +2 đến +17) với các giá trị từ -8 đến +7 cho các cặp chữ cái khác nhau. Điểm của một sắp xếp được tạo nên bởi trọng số của các cặp tương ứng trong sắp xếp trừ đi điểm phạt đối với mọi chỗ trống. Nói chung, điểm phạt chỗ trống sẽ là một hàm về chiều dài của chỗ trống. Ví dụ trên đây về số điểm bằng ma trận Dayhoff sẽ là:  Ý nghĩa đằng sau công thức tính điểm này là cách sắp xếp tối ưu hoá, điểm này có thể là cách biểu thị tốt nhất sự tương đồng về sinh học giữa hai chuỗi. Tìm ra cách sắp xếp tối ưu như vậy là nhiệm vụ chính của việc so sánh chuỗi. Thuật toán sắp xếp loại I dựa trên việc coi quá trình tìm sắp xếp tối ưu như là bài toán tìm đường đi trong một ma trận đồ thị mà trục tung và trục hoành sẽ biểu diễn 2 chuỗi theo thứ tự ban đầu của chúng. Tại mỗi điểm là giao nhau của 2 chữ cái trên 2 chuỗi thể hiện một cặp được ghép với nhau. Khi đó, mỗi một sắp xếp sẽ được biểu diễn như một đường đi trong ma trận đồ thị nói trên. Ta có thể hình dung một cách trực quan về thuật toán sắp xếp loại I dưới đây: Hình 11: Sắp xếp loại I[23] Sắp xếp loại I Cho là bảng các chữ cái giới hạn và s=(si) với i=1…|s| (|s| : số phần tử, hay lực lượng của s) và t = (ti) i=1…|t| là các phần tử xác định bởi các chữ cái trong . Định nghĩa 1 [23]: Một tập có thứ tự của các cặp (xi,yi)i=1…N,(xi,yi) {1,...,|s|}x{1,…,|t|} với x1=1Vy1=1 và xN=|s|VyN=|t| là sắp xếp loại 1 nếu: min(xi+1-xi,yi+1-yi) =1, i=1,…N-1 (1) Định nghĩa 2 [23]: Cho , là hàm khoảng trống, và là một hàm trọng số trên một cặp chữ cái trong 2 chuỗi sinh học cần so sánh, điểm số s của sắp xếp loại 1 được định nghĩa là: (2) với: (3) Nếu trong một trường hợp sắp xếp mà không bắt đầu với cặp (1,1) và không kết thúc với cặp (|s|,|t|) thì sẽ có các khoảng trống bên trái, hoặc phải của sắp xếp, chúng được gọi chung là các khoảng trống đầu cuối. Trong sắp xếp thông thường người ta không tính đến việc trừ điểm số đối với các khoảng trống đầu cuối và coi một trong 2 chuỗi như là được so sánh tương ứng với một phần của chuỗi còn lại. Nếu ta xét đến việc trừ điểm số đối với các khoảng trống đầu cuối thì hàm điểm số sẽ có dạng như sau: (4) Một hàm phạt đối với khoảng trống sẽ có dạng g(n)=a+b*(n-2) với a>b. Còn hàm trọng số thông thường đối chiếu với giá trị cho trong một ma trận trọng số, điển hình là ma trận trọng số DayHoff. Bài toán sắp xếp 2 chuỗi sinh học có thể chuyển về tìm một sắp xếp mà cho ta điểm số tối đa trong tất cả các cách sắp xếp. Thuật toán sắp xếp loại I được mô tả như một đường đi trong đồ thị định hướng loại 1 trong đó các điểm của đồ thị là tích hữu hướng {1,…,|s|} x {1,…,|t|} U {(0,0), (|s|+1,||t|+1)}. Hai điểm sau được coi là nguồn xuất phát và đích cần phải đến. Từ điểm (0,0) có thể đi tới bất kì điểm (i,j) nào thỏa mãn i=1 V j=1 và từ bất kì điểm nào có dạng (|s|,j) hay (i,|t|) đều có thể đi đến (|s|+1,|t| + 1). Nói chung mọi đường đi trong đồ thị phải thoả mãn điều kiện là số bước đi ở một trong hai hướng phải =1 hay là từ điểm (i,j) chỉ có thể đi đến được (i+1,j+k) hoặc (i+k,j+1) với k>1 và i+k < |s|+1, j+k < |t|+1. Trừ trường hợp các đường đi từ nguồn hay đến đích, tất cả các đường đi còn lại trong đồ thị e((xi,yi),(xj,yj)) sẽ được gán một trọng số như sau: (5) Trọng số cho các đường đi từ điểm nguồn hay đến điểm đích tùy thuộc vào việc ta có xét đến trừ điểm số đối với các khoảng trống bắt đầu và kết thúc hay không. Nếu có xét đến, một đường đi từ (0,0) đến (i,j) sẽ được gán trọng số là w(i,j) - G((0,0),(i,j)) và mọi đường đi đến (|s+1|,|t+1|) sẽ được gán một trọng số là –G((i,j),(|s+1|,|t+1|)). Theo nguyên tắc trên sẽ có một sự tương ứng 1-1 giữa sắp xếp loại I và đường đi trong đồ thị định hướng loại 1. Một đường đi trong đồ thị định hướng được tạo thành bởi các điểm tương ứng với một cặp chữ trên hai chuỗi sinh học. Điểm số của đồ đường đi trong đồ thị chính là điểm số của trường hợp sắp xếp trên. Bời vậy, một thuật toán chuẩn để tìm một đường đi tối ưu trong đồ thị định hướng không đối xứng với trọng số được gắn vào các cạnh khi áp dụng vào đồ thị loại 1 cũng cho ta một sắp xếp loại I tối ưu. Chúng ta sẽ xét thuật toán để tìm đường đi tối ưu này trong trường hợp không xét đến trừ điểm số đối với các khoảng trống đầu cuối. Thuật toán A-I: Trọng số của sắp xếp loại I giữa chuỗi s và t có thể được tính bằng cách lấp đầy một mảng hai chiều , L(1,j) = w(1,j) và L(i,1) = w(i,1) theo công thức: (6) Ban đầu ta sẽ lấp đầy các dòng đầu dựa và các dữ liệu có sẵn, sau đó đến dòng thứ hai dựa vào dữ liệu đã có ở dòng 1….Điểm có trọng số lớn nhất ở dòng cuối và cột cuối chính là điểm số của sắp xếp. Cách sắp xếp tương ứng với điểm số đó có thể được chỉ ra một cách chi tiết bằng phương pháp quay lui từ điểm này tới điểm mà từ đó tới nó sẽ cho giá trị lớn nhất, làm như vậy cho đến khi trở về điểm đích đầu tiên.  Chứng minh: Đỉnh vi của đồ thị được đánh số sao cho mỗi cạnh luôn trỏ tới đỉnh có chỉ số cao hơn. Độ dài lij dùng để chỉ cạnh từ vi tới vj. Và l(vi) là độ dài đường đi tối ưu từ v1 tới vi. Và: l(vi) = max (l(vj) + lij), (7) Áp dụng cách tính dựa vào bước trước với đỉnh vj. Cho các đỉnh theo thứ tự (1,1), (1,2), ..., (1, |t|), (2,1),..., (|s|, |t|). Các đỉnh trước của (i,j) là (i-k,j-1), (i-1,j-1) và (i-1,j-k). Cho đỉnh (i,j) với i=1 v j=1 đặt l(i,j) = w(sxi,tyj) và i>1^j>1 l(i - k, j - 1) + e((i - k,j - 1),(i, j)), với 2 ≤ k ≤ i – 1 l(i, j) = max l(i - 1, j - 1) + e((i - 1,j - 1),(i, j)) l(i - 1, j - k) + e((i - 1,j - k),(i, j)), với 2 ≤ k ≤ j – 1 l(i - k, j - 1) + w(si, sj) - g(k-1), với 2 ≤ k ≤ i – 1 = max l(i - 1, j - 1) + w(si, sj) l(i - 1, j - k) + w(si, sj) - g(k-1), với 2 ≤ k ≤ j – 1 l(i - k, j - 1) - g(k-1), với 2 ≤ k ≤ i – 1 = w(si, sj) + max l(i - 1, j - 1) l(i - 1, j - k) - g(k-1), với 2 ≤ k ≤ j – 1 Thuật toán trên sử dụng phương pháp qui hoạch động để lấp đầy mảng biểu diễn đồ thị. Tại mỗi điểm sau khi được lấp đầy sẽ chứa giá trị tối ưu mà từ điểm xuất phát đến nó. Như vậy, điểm cuối cùng sẽ lưu điểm số đường đi tối ưu, bằng cách áp dụng kỹ thuật quay lui ta sẽ tim được lần lượt các đỉnh mà đồ thị đi qua. Thuật toán này tương tự với thuật toán mà Needleman và Wunsch (Needleman, S.B and Wunsch, C.D J.Mol.Biol, 1970) đã đưa ra trước đó như là những người đầu tiên áp dụng phương pháp qui hoạch động vào bài toán sắp xếp chuỗi. Tại mỗi bước tiến, chương trình sẽ duyệt qua mỗi đỉnh một lần. Tại mỗi bước này, có xấp xỉ i+j giá trị được kiểm tra, trong đó i+j < |s|+|t|. Bởi vậy, khi n là số nguyên biểu diễn chiều dài của chuỗi lớn hơn trong hai chuỗi sinh học, các phép toán sẽ được thực hiện với độ phức tạp O(n3). Quá trình quay lui sẽ có độ phức tạp |s|+|t|. Chương 3: Phương pháp sắp xếp 2 chuỗi sinh học dùng Maximal unique match(MUM) Giới thiệu về MUM Mặc dù cặp gen được bảo tồn hiếm khi cùng chứa toàn bộ một chuỗi, chúng thường có chung nhiều xâu con ngắn và một số trong đó là duy nhất trong cặp gen. Xét ví dụ dưới đây, Cho 2 chuỗi S và T: S = at ga ctc a gctac t ggtcagctatt acttaccgc # T = at gt ctc t gctac ggtcagctatt c acttaccgc $ Dễ thấy 2 chuỗi S và T có chung khá nhiều xâu con: at, ctc, gctac, ggtcagctatt và acttaccgc. Trong 5 xâu con chung đó chỉ có ac không phải là duy nhất. Nó xuất hiện hơn một làn ở cả hai chuỗi. Bạn có thể chỉ ra rằng a, c, t và g cũng là xâu con chung của S và T. Tuy nhiên chúng không phải là lớn nhất vì chúng có trong ít nhất một xâu con chung dài hơn. Chúng ta chỉ xét tới những xâu con chung có độ dài lớn nhất Nhiệm vụ của chúng ta là tìm ra tất cả các xâu con chung duy nhất và có độ dài lớn nhất của hai cặp gen A và B cho trước. Mỗi xâu con chung này còn được gọi là MUM. Với hầu hết các cặp gen được bảo tồn tồn tại ít nhất một MUMs duy nhất. Vì vậy MUM có hai đặc tính: Nó xuất hiện duy nhất một lần trong cả hai cặp gen Nó không nằm trong bất kì một MUMs dài hơn (hay nó có chiều dài tối đa). Chúng ta có định nghĩa về MUM như sau: Định nghĩa 2 [22]: MUM là một xâu con chung giữa 2 gen có độ dài lớn hơn một độ dài d xác định và nó có độ dài lớn nhất, nghĩa là không thể mở rộng nó bằng cách thêm vào bất kí kí tự nào và cuối cùng MUM là duy nhất trong cả hai xâu. Với ví dụ trên giả sử d = 3, Chuỗi S và T có 4 MUMs: ctc, gctac, ggtcagctatt, acttaccgc. Xâu con ac không phải là MUM vì có độ dài bé hơn d và nó không phải là duy nhất trong cả hai chuỗi. Ta xét một ví dụ khác S = a cg a t # T = cg t a $ Cho d = 1, chúng ta sẽ có 2 MUMs: cg và t. a không phải là MUM bởi vì nó xuất hiện 2 lần ở xâu S. Một số phương pháp tìm MUM Phương pháp Brute-force: Ý tưởng của phương pháp này là trước tiên ta sẽ tìm tất cả các xâu con chung của hai chuỗi, sau đó với mỗi xâu con sẽ so sánh độ dài của nó với độ dài d cho trước và cuối cùng là xem chúng có phải là duy nhất trong cả hai chuỗi không, thuật toán được trình bày dưới đây[22]: Input: cho hai chuỗi gen S[1…m1] và T[1…m2] for every position i in S for every position i in S Find the longest common prefix P of S[i…m1] and T[j…m2] If |P| >= d and P is unique in both genomes P is MUM Output: Các MUM của hai chuỗi gen Phương pháp này có độ phức tạp về thời gian là O(m1m2), khá chậm. Phương pháp k-mers: Một phương pháp khác để tìm MUM giữa hai chuỗi x và y là tìm ra k-mers chung, k-mer đơn giản là một xâu con với độ dài k. Dưới đây là thuật toán[22]: Kí hiệu một k-mer của x là (w; 0; i) với w = xi…xi+k-1 Kí hiệu một k-mer của y là (z; 1; j) với z = yj…yj+k-1 Sắp xếp chúng theo thứ tự từ điển Đưa ra tất cả các MUM có độ dài k giữa x và y Xét cặp x, y dưới đây với k=3: x = cabbc y = bbcab 3-mers của x: (cab, 0, 1), (abb, 0, 2), (bbc, 0, 3) 3-mers của y: (bbc, 1, 1), (bca, 1, 2), (cab, 1, 3) Xắp xếp chúng: (abb, 0, 2), (bbc, 0, 3), (bbc, 1, 1), (bca, 1, 2), (cab, 0, 1), (cab, 1, 3) Đưa ra các MUM: (bbc, 0, 3), (bbc, 1, 1) và (cab, 0, 1), (cab, 1, 3) Với L là số matches. Trong trường hợp xấu nhất thời gian chạy thuật toán có độ phức tạp là O(mn) khi ta có O(m) k-mers của x và O(n) k-mers của y. Điểm yếu của thuật toán này là các MUM tìm được chưa hẳn đã là tối đa, vì thế ta lại phải tốn thêm thời gian O(L) để chọn lọc các MUM thực sự. Một giải pháp được đặt ra là tăng k và giảm L và vì thế thời gian chạy được giảm đi, tuy nhiên ta lại bỏ sót một số MUM (bé hơn k nhưng vẫn có ý nghĩa). Bởi vậy, khi đặt một giá trị cho k, chúng ta đang tạo ra một sự phù hợp giữa các MUM dài và ngắn để tối ưu thời gian chạy của thuật toán . Phương pháp sử dụng cây hậu tố: Giới thiệu về cây hậu tố Để đưa ra các MUM giữa 2 chuỗi x và y người ta sử dụng một cấu trúc dữ liệu đặc biệt là cây hậu tố. Một cây hậu tố cho xâu s chứa tất cả các suffix của s và hỗ trợ tìm kiếm với tốc độ cao, cụ thể để tìm xem 1 xâu có độ dài l có phải là xâu con của s không nó chỉ chạy trong thời gian O(l). Trước hết chúng ta phải làm rõ hai khái niệm, suffix trie và cây hậu tố. Để có một cái nhìn trực quan về suffix trie và cây hậu tố chúng ta hãy xem xét ví dụ dưới đây: Xét xâu S = abaab. a b b a a a a a b b b a baab ab Trie(abaab) Tree(abaab) Hình 12: suffix trie và cây hậu tố của xâu T = abaab Cấu trúc dữ liệu cây hậu tố, Tree(T), là một cây biểu diễn tất cả các suffix của xâu T. CTDL này có kích thước tuyến tính: |Tree(T)| = O(|T|) và có thể được xây dựng trong thời gian tuyến tính O(|T|) trong khi đó suffix trie Trie(T) lại có kích thước khá lớn|Trie(T)| = O(|T|2) bởi lẽ Trie(T) cũng biểu diễn tất cả các suffix của T tuy nhiên mỗi cạnh lại được gán nhãn là 1 kí tự, điều này làm cho kích thước của nó tăng lên khá nhiều. Định nghĩa 3[19]: Cây hậu tố T của xâu kí tự s có độ dài m là một cây định hướng có gốc mà: Cây có m lá được đánh số từ 1 tới m Mỗi đỉnh trong kể cả gốc có ít nhất 2 con Mỗi cạnh được gán nhãn bằng 1 xâu con của s Không có hai cạnh đi ra từ 1 đỉnh nào lại có nhãn cùng bắt đầu bằng một kí tự Với một lá i bất kì sự kết hợp các nhãn của các cạnh trên đường đi từ gốc tới i được gọi là 1 suffix suffix của s bắt đầu tại vị trí I xi…xm Dưới đây là 1 cây hậu tố của xâu s = xabxac: Hình 13: cây hậu tố của xâu s = xabxac Xây dựng cây hậu tố Năm 1973 Weiner đã đưa ra thuật toán cài đặt cây hậu tố với độ phức tạp tuyến tính về thời gian, cho tới năm 1976 McCreight đã phát triển một thuật toán hiệu quả hơn về không gian và năm 1995 Ukkonen đã đưa ra một thuật toán đơn giản,dễ hiểu và tối ưu hơn hai thuật toán kia, trong chương này chúng ta chủ yếu tìm hiểu thuật toán xây dựng cây hậu tố của Ukkonen. Để giảm thời gian xây dựng Cây hậu tố, Ukkonen đã đưa ra khái niệm Suffix link Định nghĩa 4 [15]: Cho trước: Xâu xα, x là một kí tự, α là một xâu Đỉnh trong v có nhãn là xα Nếu tồn tại một đỉnh s(v) với nhãn α, thì con trỏ từ v → s(v) được gọi là suffix link. Nếu α là xâu rỗng thì s(v) chính là gốc. Không có Suffix link nào xuất phát từ gốc.Trong bất kì Cây hậu tố Ti nào, nếu tồn tại đỉnh v có nhãn xα thì sẽ tồn tại đỉnh w thuộc Ti có nhãn α. Xét xâu T = t1t2 … tn$ , Pi = t1t2 … ti là prefix của T Ý tưởng của thuật toán Ukkonen: Cập nhật Trie(Pi) để thu được Trie(Pi+1) Ta xét ví dụ xây dựng Suffix Trie từ xâu abaab. a a b b a b b a a Thêm a vào Thêm b vào Thêm a vào a b b a a a a a Thêm a vào Kí tự tiếp theo được thêm vào là b a b b a a a a a Kí tự b đã xuất hiện ở đây a b b a a a a a b b b Thêm b vào cây Suffix link Hình 14: các bước tạo cây hậu tố của xâu s = abaab Để thu được Trie(Pi+1) từ Trie(Pi) : Ta sử dụng quy tắc thêm ai vào cây đã chứa ai được mô tả một cách trực quan như hình dưới Sau khi cập nhật Tại đây kí tự ai đã tồn tại => dừng cập nhật tại đây ai ai ai ai ai ai Đỉnh mới suffix links mới Suffix link Trước khi cập nhật Hình 15: Quy tắc thêm kí tự ai vào cây đã chứa ai Ý tưởng của thuật toán Ukkonen trong xây dựng cây hậu tố: construct cây hậu tố T1 For i = 1 to m – 1 { phase i+1 – buid Ti+1 from Ti by adding character number i+1} Có 3 quy tắc khi thêm một kí tự vào cây cho trước: Áp dụng với cây hậu tố tạo bởi xâu S = axabx 4 a a x b a x x x 2 1 b x b b x 3 Hình 16: Cây hậu tố T của xâu S = axabx Tại bước i+1, extention thứ j đặt suffix S[j…i] là β Nếu β kết thúc tại lá trong cây Ti → thêm S(i+1) vào cuối nhãn của cạnh ứng với lá Áp dụng vào cây T ở hình 5: Thêm vào kí tự b, ta được xâu S = axabxb 4 a a x b a x x x 1 b x b b x 3 b b b b 2 Hình 17: Cây hậu tố T của xâu S = axabxb Nếu không có đường nào từ kết thúc của β bắt đầu với S(i+1)nhưng có ít nhất một đường được gán nhãn tiếp tục từ kết thúc của β→ Tạo cạnh ứng với lá mới bắt đầu từ kết thúc của β gán cho cạnh đó nhãn S(i+1). Tạo ra đỉnh mới nếu β kết thúc bên trong 1 cạnh, gán cho lá ứng với điểm cuối của cạnh này là j. 4 a a x b a x x x 2 1 b x b b x 3 b b b b b 5 Trở lại ví dụ ở hình 6: Điểm cuối của β là x; không có đường nào chứa “xb” nhưng có 1 đường được gán nhãn chứa “xa” → extention j=5 được tạo. Hình 18: Cây hậu tố T của xâu S = axabxb Nếu tồn tại đường đi từ kết thúc của β mà bắt đầu với S(i+1)→ không cần làm gì nữa. Thủ tục xây dựng Suffix Trie[7] Input: Xâu T = t1t2 … tn$ Output: Trie(T) Kí hiệu: son(v,α) = w Nếu có một arc từ v tới w có nhãn α, chú ý: son(v,ε) = v Construct Trie(t1): include root, node v, arc son(root, t1) = v, suffix links slink(v) := root, slink(root) := root for i := 2 to n do begin vi-1 := leaf of Trie(t1…ti-1) of string t1…ti-1 (the deepest leaf) v := vi-1; v´ := 0 while node v haven’t any arc from ti do begin creat new node v´´ and new arc son(v,ti) = v´´ if v´ ≠ 0 then slink(v) := v´´ v := slink(v); v´ := v´´ end for node v´´ and v´´= son(v,ti) do if v´´ = v´ then slink(v’) := root else slink(v´) := v´´ Thủ tục xây dựng Cây hậu tố[7] Input: Xâu T = t1t2 … tn$ Output: Tree(T) Kí hiệu: son(v,α) = w Nếu có một arc từ v tới w có nhãn α, chú ý: son(v,ε) = v Ukkonen đưa ra Hàm Canonize(v, α) [7]: while son(v, α´) ≠ 0 : α = α´ α´´, | α´| > 0 do v := son(v, α´); α := α´´ return (v, α) 1. Create Tree(t1); suffix link slink(root) := root (v, α) := (root, ε) /* (v, α) là đỉnh bắt đầu */ 2. for i := 2 to n+1 do 3. v´ := 0 4. while not exit any arc which start from v has prefix label αti do 5. if α ≠ ε then /* chia arc w = son(v, αη) thành hai */ 6. son(v, α) := v´´; son(v´´,ti) := v´´´; son(v´´,η) := w 7. else 8. son(v,ti) := v´´´; v´´ := v 9. if v´ ≠ 0 then slink(v´) := v´´ 10. v´ := v´´; v := slink(v); (v, α) := Canonize(v, α) 11. if v´ ≠ 0 then slink(v´) := v 12. (v, α) := Canonize(v, αti) /* (v, α) = đỉnh tiếp theo của vòng lặp kế */ Thuật toán xây dựng cây hậu tố của Ukkonen có độ phức tạp về thời gian là O(m) với m là độ dài của xâu nguồn, ta có thể kiểm tra xem một xâu có độ dài n có phải là xâu con của xâu nguồn không trong thời gian O(n). Phương pháp tìm MUM sử dụng cây hậu tố Bài toán: Đưa ra tất cả các MUM giữ hai xâu S1, S2 cho trước. Input: Hai xâu S1 và S2; Output: Tất cả các MUM của hai xâu; Ý tưởng cơ bản của thuật toán là xây dựng cây hậu tố cho xâu có độ dài bé hơn trong hai xâu nhằm giảm thời gian duyệt cây do đó tăng tốc độ tìm kiếm các MUM. Sau đó duyệt xâu còn lại trên cây vừa xây dựng được để đưa ra được các xâu con chung, cuối cùng ta mở rộng các xâu con này ra tới mức tối đa và thu được các MUM. Giả sử |S1|<|S2|, ta sẽ xây dựng cây hậu tố T cho xâu S1. Quy ước: Cấu trúc MUM lưu trữ các thông tin sau (start1, start2, length). Start1 là vị trí bắt đầu của MUM ở xâu S1 và start2 là vị trí bắt đầu của MUM ở xâu S2, length là độ dài của MUM, length >= k. k là độ dài tối thiểu của 1 MUM, S1(p) là xâu con của S1 từ vị trí p đến p+k-1, S2(p) là xâu con của S2 từ p đến p-k+1, SS là mảng lưu trữ các MUM hiện tại, len1: độ dài S1, len2: dộ dài S2. 1 MUM sa = (start1, start2, length) trong danh sách SA được gọi là mở rộng được (extendable), nếu (start1 + p - start2) là 1 giá trị trong danh sách L1 Thủ tục tìm MUM sử dụng cây hậu tố: Build suffix tree T1 for S1, set SA = [], L1=[], SS=[]; for p := 0 to len2 - k do slide S2(p) on T1; if S2(p) end at a node of T1 then Save all leaf behind this node to list L1; /*L1 lưu tất cả các vị trí xuất hiện của S2(p) trong xâu S1*/ Save all MUM to SS; While all MUM in SS isn’t extentable if MUM sa= (start1, start2, length) is extentable then length++; delete value (start1 +1) in L1; else save sa to SA; if L1!=NULL then save all value in L1 to SA; Để đánh giá khoảng cách của hai chuỗi sinh học, chúng tôi đề xuất ra một tiêu chuẩn dựa trên số lượng và độ dài các MUM vừa tìm được: Giả sử với hai chuỗi cho trước S1, S2 như trên Chúng ta tìm được: m1 MUM có độ dài k m2 MUM có độ dài k+1 …. mi MUM có độ dài k+j với j <= MIN(len1, len2) Khi đó khoảng cách giữa hai chuỗi S1, S2 sẽ được xác định bởi công thức: (8) Khi D càng bé tức là 2 chuỗi S1, S2 có càng nhiều MUM thì khoảng cách giữa chúng càng bé. Chương 4: Đánh giá kết quả thực nghiệm Kết quả thực nghiệm Song song với phương pháp sử dụng cây hậu tố để sắp xếp 2 chuỗi gen, chúng tôi có thử nghiệm trên thuật toán quy hoạch động truyền thống với cùng một bộ dữ liệu là các chuỗi gen của chủng virus truyền bệnh đậu mùa (poxvirus). Cụ thể là 12 chuỗi gen trong file pox_subsets.fas từ trang web Bảng 6: Dữ liệu sử dụng trong thực nghiệm STT Tên virus Độ dài 1 CMLV-CMS Camelpox virus CMS 202205 2 CMLV-M96 Camelpox virus isolate M-96 from Kazakhstan 205719 3 Cowpox virus strain Brighton Red 224499 4 CPXV-GRI Cowpox virus strain GRI-90 223666 5 ECTV-MOS Ectromelia virus strain Moscow 209771 6 ECTV-NAV Ectromelia virus strain Nava 207620 7 GTPV-G20LKV Goatpox virus strain G20-LKV 149723 8 GTPV-Pellor Goatpox virus strain Pellor 149599 9 LSDV-1959 Lumpy skin disease virus isolate Neethling vaccine LW 1959 150509 10 LSDV-NEE Lumpy skin disease virus strain Neethling isolate 2490 150733 11 LSDV-WARM Lumpy skin disease virus isolate Neethling Warmbaths LW 150793 12 MPXV-ZRE Monkeypox virus strain Zaire-96-I-16 196858 Sau khi tính khoảng cách của từng cặp trong 12 chuỗi gen theo phương pháp quy hoạch động, chúng tôi thu được kết quả như sau. Bảng 7: Khoảng cách của các cặp chuỗi theo quy hoạch động STT 1 2 3 4 5 6 7 8 9 10 11 1 2 3694 3 37164 34751 4 33928 33167 14603 5 30929 30387 29754 28535 6 28654 28577 29431 28132 4603 7 91236 93969 107616 106953 96544 95070 8 91275 94005 107668 107006 96610 95107 130 9 91059 93801 107397 106746 96355 94916 4978 5086 10 90975 93719 107318 106623 96249 94828 4533 4647 2386 11 90973 93709 107314 106622 96257 94840 4551 4667 2412 251 12 27076 28934 38489 34523 32944 31311 87338 87383 87162 87051 87048 Với các giá trị của k (độ dài tối thiểu của MUM) lần lượt là 10, 15, 20 chúng tôi thu được các kết quả sau. Bảng 8: Khoảng cách của các cặp chuỗi theo phương pháp sử dụng cây hậu tố (k=10) STT 1 2 3 4 5 6 7 8 9 10 11 1 2 43.8042 3 1754.69 1630.38 4 1448.55 1371.24 1384.13 5 1836.83 1758.72 1685.76 1550.98 6 1837.52 1729.34 1623.37 1539.7 42.0266 7 2524.29 2486.64 2486.86 2590.26 2509.4 2467.77 8 2523.27 2485.51 2485.67 2589.05 2507.9 2466.34 6.11453 9 2580.97 2549.11 2538.85 2648.54 2569.18 2528.26 688.693 688.47 10 2590.17 2556.87 2545.11 2658.31 2574.59 2534.24 643.652 643.113 438.003 11 2591.77 2558.43 2547.32 2660.25 2575.6 2535.55 646.444 646.01 440.854 71.7801 12 1718.1 1721.51 1703.24 1338.2 1752.16 1744.39 2518.23 2517.19 2571.79 2583.63 2585.4 Bảng 9: Khoảng cách của các cặp chuỗi theo phương pháp sử dụng cây hậu tố (k=15) STT 1 2 3 4 5 6 7 8 9 10 11 1 2 44.3114 3 3299.34 2947.05 4 2314.57 2131.55 2186.09 5 3597.84 3361.81 3126.16 2618.53 6 3636.28 3332.38 3007.95 2616.43 42.5106 7 253501 244609 257640 260283 251289 237776 8 253291 244407 257427 259551 251081 237579 6.14611 9 250183 234064 253822 271364 249973 235717 1565.6 1564.45 10 259914 251533 255903 269404 251983 236620 1352.18 1349.72 667.714 11 253759 237230 255128 270037 250790 236071 1367.77 1365.7 674.962 76.0836 12 3152.46 3179.8 3140.09 2050.75 3295.6 3301.45 253309 253099 259460 264730 265371 Bảng 10: Khoảng cách của các cặp chuỗi theo phương pháp sử dụng cây hậu tố (k=20) STT 1 2 3 4 5 6 7 8 9 10 11 1 2 44.3147 3 3370.56 3067.47 4 2338.75 2152.66 2221.18 5 3671.14 3424.2 3225.93 2655.09 6 3712.42 3395.16 3091.39 2659.13 42.5147 7 1269183 1291239 1483626 1441598 1277991 1264886 8 1268132 1290170 1482398 1440404 1276932 1263839 6.1464 9 1299988 1322580 1507047 1378096 1278354 1265246 1593.66 1592.33 10 1306958 1329671 1585794 1408418 1335367 1321674 1372.93 1370.05 672.102 11 1307131 1329847 1586004 1406080 1308799 1295378 1388.9 1386.45 679.429 76.1303 12 3211.31 3242.54 3205.8 2069.62 3360.46 3367.13 1151123 1150170 1166047 1212220 1212381 Để so sánh thời gian sắp xếp giữa hai thuật toán chúng tôi sử dụng hệ thống máy tính Pentum 4 tốc độ 3Ghz với bộ nhớ trong 512MB để kiểm tra lần lượt từng cặp gen trong số 12 mẫu gen kể trên., thời gian cụ thể được trình bày ở bảng 11 Bảng 11: Thời gian chạy của hai thuật toán Phương pháp Tổng số test Thời gian chạy Thời gian chạy 1 test Quy hoạch động 66 36 giờ 33 phút 43 giây Sử dụng cây hậu tố 66 2 phút 12s 2 giây Đánh giá kết quả Để đánh giá kết quả thực nghiệm thu được từ phương pháp quy hoạch động truyền thống và phương pháp sử dụng cây hậu tố, chúng tôi đánh giá độ tương tự giữa hai ma trận xây dựng được từ cả hai phương pháp bằng cách sau: Xét một bộ 3 xâu vào a,b,c nếu ở ma trận khoảng cách xây dựng được bằng phương pháp quy hoạch động D(a,b)< D(a,c) thì ở ma trận khoảng cách xây dựng được bằng phương pháp sử dụng cây hậu tố D(a,b)< D(a,c), nghĩa là thứ tự của các xâu vẫn giữ nguyên ở cả hai ma trận. Hai ma trận có càng nhiều bộ ba như thế này thì độ tương tự càng cao. Trở lại kết quả thực nghiệm trên, với k=10, chúng tôi tính ra được độ tương tự giữa hai ma trận là 80,43%. Với k=15, độ tương tự giữa hai ma trận là 84,47%. Và với k=20, độ tương tự giữa hai ma trận là 86,49%. Nghĩa là với k = 20 thì ma trận xây dựng được từ phương pháp sử dụng cây hậu tố giống nhất với ma trận xây dựng bằng phương pháp quy hoạch động Trong khi đó phương pháp sắp xếp sử dụng cây hậu tố chạy nhanh hơn phương pháp quy hoạch động truyền thống hơn 1000 lần, đây là một cải thiện về tốc độ rất đáng kể. Kết luận Để tìm ra mối quan hệ tiến hóa giữa 2 loài sinh vật, người ta nghiên cứu về sắp xếp các chuỗi gen của 2 loài đó. Quá trình nghiên cứu này đòi hỏi vận dụng một số kiến thức về sinh học phân tử và các thuật toán sắp xếp. Sau đó, thông qua thực nghiệm, vấn đề sẽ được làm rõ. Trong khóa luận này, chúng tôi đã sử dụng cây hậu tố để tìm ra các MUM giữa 2 chuỗi gen, từ đó, xây dựng ma trận khoảng cách giữa chúng. Kết quả thực nghiệm cho thấy, ma trận khoảng cách xây dựng từ phương pháp sử dụng cây hậu tố và ma trận khoảng cách có được từ phương pháp quy hoạch động có độ tương đồng khá cao. Đặc biệt, tốc độ nhanh hơn 1000 lần, đây là một cải thiện rất đáng kể về tốc độ so với phương pháp quy hoạch động. Trong tương lai, chúng tôi sẽ nghiên cứu đưa ra một tiêu chuẩn mới để so sánh khoảng cách giữa 2 chuỗi dựa trên các MUM tìm được để có được kết quả sắp xếp tối ưu nhất. Sau đó, sẽ xây dựng cây phát sinh loài nhằm chỉ ra mối quan hệ tiến hóa giữa các loài sinh vật, giải quyết một trong những bài toán lớn của Tin sinh học hiện nay. Tài liệu tham khảo Hồ Tú Bảo, Phạm Thọ Hoàn, School of Knowledge Science Japan Advanced Institute of Science and technology, Tin sinh học Khái niệm và bài toán cơ bản Một vài kết quả nghiên cứu, 2005, tr. 1-13. Lê Đình Lương, Phan Cự Nhân, Cơ sở di truyền học, NXBGD, 1997. Aluízio Borém, Fabrício R. Santos, David E. Bowen, Understanding Biotechnology, Prentice Hall, 2003. Damian Counsell, Bioinformatics FAQ, 2003. Divya R. Singh, Faster Sequence Alignment using suffix tree and Data-Mining Techniques, February 2005. Esko Ukkonen, University Helsinki, suffix tree and suffix array techniques for pattern analysis in strings, Erice School 30 Oct 2005. Esko Ukkonen, On-line construction of suffix tree, 1995. Esti Yeger Lotem, Gideon Greenspan, Lecture I: Introduction & Text Based Search, 2001. Geoffrey Alan Mazeroff, Probabilistic Suffix Models for Windows Application Behavior Profiling: Framework and Initial Results, December 2004. Hong T. Lam , Computational Genomics: Sequence Alignment , January 6, 2004. Marc Rehmsmeier, database searching with phylogenetic trees, October 2001. Mark S. Boguski, Trends Guide to Bioinformatics, 1998. Max Planck Institute for Molecular Genetics, Germany, Online Lectures on Bioinformatics. Misha Taylor, developing a bioinformatics utility belt to eliminate search redundancy from the ever-growing databases, spring semester 2003. Nathan Edwards, Lecture 12: suffix tree, Algorithms in Biosequence Analysis - Fall, 2005. Peter Clote and Rolf Backofen, Ludwig-Maximilians-Universit¨at M¨unchen, Computational molecular biology: an introduction,JOHN WILEY & SONS, tr.1-22. Pierre Baldi, Soren Brunak, Second edition of Bioinformatics, The machine Learning approach, Masachusettles institute of technology, 2001. S. Schulze-Kremer, Advances in Molecular Bioinformatics, IOS press 1994. Saad Mneimneh, Computational biology lecture 18, Stefan Kurtz, Lecture notes for a course in the Winter Semester 2000/2001 Foundations of Sequence Analysis, July 18 2002. T K Attwood and D J Parry-Smith, Introduction to Bioinformatics, Prentice-Hall 1999 [Longman Higher Education; ISBN 0582327881]. Wing-Kin Sung Scribe: Hady Gunawan, Steven Halim, Tan Soon Heng, Combinatorial Methods in Bioinformatics 2004/2005 Semester 1, Lecture 3: Whole Genome Alignment - August 27, 2004. Online Lectures on Bioinformatics.