TÓM TẮT LUẬN VĂN
Cây sinh loài mô tả lịch sử tiến hóa của một nhóm các loài 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 tính của mỗi loài được chúng ta quan tâm ở đây tương
ứng với các bộ gen. Gen là các chuỗi DNA được bao gồm từ các kí tự A, G, C và T
hợp thành. Cây sinh loài là một cây mà các nút lá (taxa) của nó có thể là các vật sống
hiện tại ngày nay, các nút trong của cây đó là các tổ tiên của các nút lá. Tái cấu trúc
cây sinh loài chính là tìm những gen phù hợp nhất để đưa vào các nút tổ tiên hoặc là
đưa ra một cây sinh loài phù hợp nhất để giải thích quá trình tiến hoá.
Tuy nhiên, việc nghiên cứu cây sinh loài cho nhiều hướng tiếp cận. Mỗi phương
pháp có những ưu điểm và khuyết điểm của nó. Phương pháp ước lượng hợp lý cực
đại được chọn ở đây là phương pháp phức tạp nhất nhưng lại là phương pháp cho kết
quả tin cậy nhất. Công cụ chính sử dụng trong phương pháp này là Đại số thống kê và
Đại số máy tính. Đó là những lãnh vực phát triển mạnh mẽ trong những năm gần đây.
Thống kê là ngành khoa học phân tích dữ liệu. Đối với các chuỗi DNA thì
thống kê sẽ xây dựng những mô hình quá trình phát sinh dữ liệu. Đưa ra những kết
luận chung về quá trình phát sinh đó. Mô hình thống kê là nguyên tắc cơ bản đối với
các gen. Đại số thống kê làm sáng tỏ cho những ý tưởng trọng tâm về phân tích dữ liệu
rời rạc nói riêng và phân tích chuỗi sinh học nói riêng.
Ước lượng hợp lý cực đại (Maximum Likelihood Estimation – MLE) được
công thức hoá trong Xác suất cổ điển, nó có tính chất của một ước lượng tốt. Phương
pháp MLE đánh giá những tham số của một mô hình thối lui. MLE dẫn đến việc giải
quyết là làm cực đại tích của những đa thức.
Đại số máy tính là một lãnh vực mới, nó cung cấp những nền tảng để giải bài
toán MLE trên máy tính.
Đề tài này tập trung vào việc nghiên cứu mô hình xác suất thống kê trên cây
sinh loài từ những dữ liệu là các gen của sinh vật sống. Sau đó sử dụng những nền tảng
toán học, đại số máy tính để giải quyết bài toán hợp lý cực đại của mô hình xác suất
trên. Mục tiêu cuối cùng là tìm một cây sinh loài thích hợp nhất để giải thích sự tiến
hoá. Những kết quả của luận văn đã làm như sau:
- Về phương pháp: Chọn phương pháp đáng tin cậy nhất là phương pháp ước
lượng hợp lý cực đại cho mô hình hóa bài toán. Giải phương trình hợp lý bằng
phương pháp tính toán đại số để tìm kết quả chính xác.
- Về tính toán: Viết một chương trình để mô hình hóa ước lượng hợp lý cực đại
trên cây sinh loài và chạy tìm nghiệm phương trình hợp lý trên một số cây sinh
loài nhỏ 3 và 4 taxa ở một số mô hình.
MỤC LỤC
LỜI CAM ĐOAN 1
LỜI CẢM ƠN 2
TÓM TẮT LUẬN VĂN 3
DANH MỤC BẢNG 4
DANH MỤC HÌNH .5
MỤC LỤC 6
Chương 1. GIỚI THIỆU ĐỀ TÀI 9
1.1. Giới thiệu 9
1.2. Cấu trúc luận văn 10
Chương 2. CƠ SỞ LÝ THUYẾT VỀ CÁC CẤU TRÚC ĐẠI SỐ VÀ XÁC SUẤT
THỐNG KÊ 12
2.1. Một số cấu trúc đại số cơ bàn .12
2.1.1. Lý thuyết nhóm .12
2.1.2. Lý thuyết vành 13
2.1.3. Trường 14
2.1.4. Vành đa thức .14
2.1.5. Ma trận 15
2.1.6. Định thức 15
2.1.7. Không gian vector .16
2.1.8. Đa tạp đại số .18
2.2. Các khái niệm về xác suất thống kê .18
2.2.1. Định nghĩa về xác suất 18
2.2.2. Xác suất có điều kiện 19
2.2.3. Đại lượng ngẫu nhiên và hàm phân phối 20
2.2.4. Các đặc trưng của đại lượng ngẫu nhiên .20
2.2.5. Lý thuyết mẫu .21
2.2.6. Ước lượng tham số 22
2.2.7. Sơ lược về ước lượng hợp lý cực đại 22
Chương 3. ƯỚC LƯỢNG HỢP LÝ CỰC ĐẠI TRÊN MẪU QUAN SÁT 25
3.1. Ước lượng hợp lý cực đại là gì? 25
3.1.1. Đặt vấn đề .25
3.1.2. Khái quát về ước lượng hợp lý cực đại .25
3.1.3. Ví dụ về ước lượng hợp lý cực đại .26
3.2. Giải bài toán ước lượng hợp lý cực đại 26
3.2.1. Nguyên lý ước lượng hợp lý cực đại 26
3.2.2. Logarit hàm hợp lý 26
3.3. Tổng quát hóa bài toán ước lượng hợp lý cực đại 27
3.3.1. Ước lượng hợp lý cực đại trên mẫu quan sát 27
3.3.2. Một số phương pháp giải phương trình hợp lý .28
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 .30
4.1. Giới thiệu sơ lược về cây sinh loài 30
4.2. Các nghiên cứu phát sinh sinh loài .31
4.3. Mô hình ước lượng hợp lý cực đại trên cây sinh loài .32
4.4. Mô hình tiến hóa 33
Chương 5. BẤT BIẾN TRÊN CÂY SINH LOÀI .37
5.1. Dẫn nhập .37
5.2. Mô hình xác suất trên cây sinh loài 38
5.2.1. Mô hình bài toán cây sinh loài 38
5.2.2. Nhóm Abel và sự liên hệ với các ma trận chuyển đổi 39
5.3. Biến đổi Fourier .40
5.4. Toạ độ Fourier .42
5.5. Áp dụng tìm bất biến trên một cây sinh loài .42
5.5.1. Mô hình bài toán .42
5.5.2. Các khả năng xảy ra trên các nút lá 43
5.5.3. Các lớp xác suất tương đương 43
5.5.4. Chuyển đổi Fourier .44
5.5.5. Kết quả tìm được 45
5.6. Những tính chất của thành phần bất biến .46
Chương 6. GIẢI PHƯƠNG TRÌNH HỢP LÝ 47
6.1. Quỹ tích hợp lý trên một đa tạp .47
6.2. Ma trận Jacobi của các đa thức bất biến .47
6.2.1. Gradient- Vector vận tốc .47
6.2.2. Ma trận Jacobi của các đa thức bất biến .48
6.2.3. Không gian tiếp xúc 49
6.3. Bài toán cực trị điều kiện 49
6.4. Bậc của hợp lý cực đại .50
6.5. Các thuật toán 50
6.6. Áp dụng giải phương trình hợp lý 51
Chương 7. CHƯƠNG TRÌNH THỰC HIỆN .53
7.1. Sơ đồ khối chương trình 53
7.2. Sơ lược về chương trình 54
7.3. Kết quả chương trình 54
Chương 8. TỔNG KẾT – ĐÁNH GIÁ 57
8.1. Tổng kết 57
8.2. Những đóng góp của luận văn 57
8.3. Hướng phát triển .58
TÀI LIỆU THAM KHẢO .59
76 trang |
Chia sẻ: lvcdongnoi | Lượt xem: 2623 | Lượt tải: 1
Bạn đang xem trước 20 trang tài liệu 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ỏ, để xem tài liệu hoàn chỉnh bạn click vào nút DOWNLOAD ở trên
r đối với những tổ hợp của xác suất trên đối với
nhóm . Để làm điều đó, chúng ta thay thế phân bố gốc nG :π →G R bởi một hàm số
mới : nπ → G R như sau:
1 1 2
1
( )
( ,..., )
0
neáu
ngöôïc laïi
n
n
h h h ... h
h h
ππ = = =⎧= ⎨⎩
Từ đó chúng ta có:
∑ ∏
∈ =
−=
n
n Ghh
n
i
ii
i
nn ghfhhggp
),...,( 1
)(
11
1
)(),...,(~),...,( π
vì thế p là chập của hai hàm số trên . Biến đổi Fourier cho ra: nG
∏
=
=
n
i
i
i
nn fq
1
)(
11 )(ˆ),...,(~ˆ),...,( χχχπχχ
bởi công thức tính chập là độc lập của f(i) trong biến đổi Fourier. Mặt khác
∑
∑
=〉〈=
〉〈=
∈
∈
)...(ˆ)(.),...
),...,(~.),...,(),,...,(),...,(~ˆ
2121
111),..,(1 1
nnGg
nnnGggn
gg
ggggn
n
χχχππχχχ
πχχχχπ
cho nên
∏
=
=
n
i
i
i
nn fq
1
)(
11 )(ˆ),...,(ˆ),...,( χχχπχχ
Ví dụ trên là cơ sở để giới thiệu sự cần thiết để chứng minh những kết quả tổng
quát sau.
Định lý ([Evans and Speed, 1993]): Cho là phân phối có điều
kiện của một mô hình cơ sở nhóm đối với cây sinh loài T được giới thiệu ở phần trên.
Thì biến đổi Fourier của p có dạng
),...,( 1 nggp
)(ˆ),...,(ˆ),...,(
}{\)( )(
)(
11 ∏ ∏
∈ Λ∈
=
rTVv vl
l
v
nn fq χχχπχχ
với )(vΛ là tập của lá có v như là cha.
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 42
Thay thế tọa độ gốc bởi tọa độ Fourier , kết quả của là các
đơn thức của các tham số.
niiip ..21 niiiq ..21 niiiq ..21
5.4. Toạ độ Fourier
Mỗi một tọa độ Fourier của 2n hoặc 4n tọa độ được ký hiệu bởi . Chú ý,
với phân phối tại gốc là phân phối đều và mô hình chúng ta đang xét có cấu trúc nhóm
như Jukes-Cantor hay Kimura-2, Kimura-3, biến đổi Fourier từ theo và
ngược lại theo định lý trên (ở phần 5.3) như sau:
niiiq ..21
niiip ..21 niiiq ..21
∑=
n
n
n
n
jj
jjn
jj
iii qiip
,...,
...1..
1
1
1
21
,)()...( χχ
∑=
n
n
n
n
jj
jjn
ii
niii pjjk
q
,...,
...1..
1
1
1
21
.)()...(1 χχ
Ở đây iχ là đặc trưng của nhóm kết hợp đến phần tử thứ i của nhóm. Bảng đặc
trưng của nhóm chúng ta sử dụng là và 2Z 2 2⊕Z Z sau:
111
110
1 0
-
và
1111
1111
1111
1111
- -T
- - G
- -C
A
T G C A
Nói cách khác, là chính là phần tử (i, j) tương ứng trong bảng đặc trưng
trên.
)( jiχ
5.5. Áp dụng tìm bất biến trên một cây sinh loài
5.5.1. Mô hình bài toán
Để làm rõ vấn đề trên chúng ta xem xét một cây sinh loài cụ thể sau: Cho cây
sinh loài có gốc với 3 lá với hình 8 sau:
Hình 8: Cây sinh loài có gốc với 3 nút lá
Để đơn giản, giả sử mô hình chúng ta đang xét là mô hình Jukes – Cantor với 4
trạng thái. Ngoài ra, các kí tự tại nút gốc có phân phối đều, nói cách khác xác suất xuất
hiện các ký tự {A, G, C, T} tại gốc là bằng nhau và bằng 0.25.
Theo hình dáng cây trên, chúng ta có 3 cạnh bằng nhau, cho nên ma trận xác
xuất chuyển đổi như 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 43
0 1 1 1
1 0 1 1
1 1 0 1
1 1 1 0
c c c c
c c c c
c
c c c c
c c c c
⎛ ⎞⎜ ⎟⎜ ⎟= ⎜ ⎟⎜ ⎟⎝ ⎠
(Với ma trận trên ta có ràng buộc c0+3c1=1).
5.5.2. Các khả năng xảy ra trên các nút lá
Như trên, vì cây có 3 lá, nên số trường hợp có thể xảy ra trên lá là 64 khả năng.
(Xem Phụ lục 1)
5.5.3. Các lớp xác suất tương đương
Trong 64 toạ độ trên, thì chỉ có 4 lớp tương đương sau:
Lớp 1 có 4 toạ độ:
3 3
0 1
1 3
4 4AAA GGG CCC TTT
p p p p c= = = = + c
Lớp 2 có 36 tọa độ:
2
0 1 0
1 1
4 4
AAC AAG AAT ACA ACC AGA AGG ATA ATT
CAA CAC CCA CCG CCT CGC CGG CTC CTT
GAA GAG GCC GCG GGA GGC GGT GTG GTT
TAA TAT TCC TCT TGG TGT TTA TTC TTG
p p p p p p p p p
p p p p p p p p p
p p p p p p p p p
p p p p p p p p p
c c c
= = = = = = = =
= = = = = = = = =
= = = = = = = = =
= = = = = = = = =
= + 2 31 112c c+
Lớp 3 có 24 tọa độ:
2 3
0 1 1
3 1
4 4
ACG ACT AGC AGT ATC ATG CAG CAT CGA
CGT CTA CTG GAC GAT GCA GCT GTA GT
TAC TAG TCA TCG TGA TGC
p p p p p p p p p
Cp p p p p p p p p
p p p p p p c c c
= = = = = = = =
= = = = = = = = =
= = = = = = = +
Đặt 0 1, 2,p p p là tổng các xác suất trong một lớp tương đương trên, ta có:
3 3
0 0 3p c c= + 1
3
1
1c
m
2 2
1 0 1 0 19 9 18p c c c c c= + +
2 3
2 0 118 6p c c= +
Vậy ánh xạ chúng ta cần xem xét , trong đó: : df →^ ^
Trong đó d = 2 (c0 và c1) và cây có 3 lớp đại diện cho 64 khả năng xảy ra ở các
nút lá cho nên m = 3.
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 44
Vì thế thực tế ánh xạ chúng ta đang xét đó là:
2 3:f →\ \
5.5.4. Chuyển đổi Fourier
Sử dụng công thức chuyển đổi Fourier
∑=
n
n
n
n
jj
jjn
ii
niii pjjk
q
,...,
...1..
1
1
1
21
.)()...(1 χχ
ta được:
3 2 2 3
0 1 2 0 0 1 0 1 1
3 2 2
0 1 2 0 0 1 0 1
9 27 27
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
AAA
AAG
q p p p c c c c c c
q p p p c c c c c
= + + = + + +
= + = + + 31c
3 2 2
0 1 2 0 0 1 0 1
3 2 2
0 1 2 0 0 1 0 1
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
AAC
AAT
q p p p c c c c c
q p p p c c c c c
= + = + +
= + = + +
3
1
3
1
c
c
3 2 2
0 1 2 0 0 1 0 1
3 2 2
0 1 2 0 0 1 0 1
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
AGA
AGG
q p p p c c c c c
q p p p c c c c c
= + = + +
= + = + +
3
1
3
1
c
c
3 2 2 3
0 1 2 0 0 1 0 1
3 2 2 3
0 1 2 0 0 1 0 1
1 1 1 1- -
18 18 2 2
1 1 1 1- -
18 18 2 2
1 1
6 6
1 1
6 6
AGC
AGT
q p p p c c c c c
q p p p c c c c c
= + = +
= + = +
1
1
-
-
1
6
1
6
c
c
3 2 2
0 1 2 0 0 1 0 1
3 2 2
0 1 2 0 0 1 0 1
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
1 1 1 1- -
18 18 2 2
1 1
6 6
ACA
ACG
q p p p c c c c c
q p p p c c c c c
= + = + +
= + = +
3
1
3
1-
1
6
c
c
3 2 2
0 1 2 0 0 1 0 1
3 2 2 3
0 1 2 0 0 1 0 1
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
1 1 1 1- -
18 18 2 2
1 1
6 6
ACC
ACT
q p p p c c c c c
q p p p c c c c c
= + = + +
= + = +
3
1
1-
1
6
c
c
3 2 2
0 1 2 0 0 1 0 1
3 2 2 3
0 1 2 0 0 1 0 1
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
1 1 1 1- -
18 18 2 2
1 1
6 6
ATA
ATG
q p p p c c c c c
q p p p c c c c c
= + = + +
= + = +
3
1
1-
1
6
c
c
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 45
3 2 2 3
0 1 2 0 0 1 0 1
3 2 2
0 1 2 0 0 1 0 1
1 1 1 1- -
18 18 2 2
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
1 1
6 6
ATC
ATT
q p p p c c c c c
q p p p c c c c c
= + = +
= + = + +
1
3
1
- 1
6
c
c
3
1c
Tất cả q còn lại đều mang giá trị 0.
Các q trên được chia thành 3 lớp tương đương:
Lớp 1(có một tọa độ):
3 2 2
1 2 3 0 0 1 0 19 27 27 AAAq p p p c c c c c= + + = + + +
Lớp 2 (có 9 tọa độ):
3 2 2 3
1 2 3 0 0 1 0 1 1
1 1 1 1 1 5 1- -
9 81 27 9 9 9 3
AAG AAC AAT AGA AGG ACA ACC ATA ATTq q q q q q q q q
p p p c c c c c c
= = = = = = = =
= + = + +
Lớp 3 (có 6 tọa độ):
3 2 2 3
1 2 3 0 0 1 0 1
1 1 1 1- -
18 18 2 2
1 1
6 6
AGC AGT ACG ACT ATG ATCq q q q q q
1-
1
6
p p p c c c c c
= = = = =
= + = + c
3
1
Đặt là tổng giá trị của từng lớp tương đương trên, thì 0 1 2, ,q q q
3 2 2
0 0 1 2 0 0 1 0 19 27 27 q p p p c c c c c c= + + = + + +
3 2 2
1 0 1 2 0 0 1 0 1
1 1- - 5
9 3
q p p p c c c c c c= + = + + 313
3 2 2
2 0 1 2 0 0 1 0 1
1 1- - 3 3
3 3
q p p p c c c c c c= + = + 31-
5.5.5. Kết quả tìm được
Bất biến cần tìm là được:
2 3
0 2 1- 0q q q =
Hay bất biến cần tìm:
2 3
0 1 2 0 1 2 0 1 2
2 2 2 3 2 2 3
0 2 0 1 0 2 1 1 2 1 2 2
1 1 1 1( )( - ) - ( - ) 0
3 3 9 3
8 16 4 80 8 4 4- - -
3 27 9 729 81 27 27
p p p p p p p p p
p p p p p p p p p p p p
+ + + + =
⇔ + + + 0=
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 46
5.6. Những tính chất của thành phần bất biến
Cũng theo các tác giả [Evans and Speed, 1993] thì với phép biến đổi Fourier
trên, chúng ta sẽ tìm được tất cả các thành phần bất biến trên cây sinh loài. Và một
điều quan trọng nữa là các thành phần bất biến trên là những đa thức thuần nhất.
Thành phần bất biến tầm thường nhất là 1i
i
p =∑ mà ta đã biết. Những thành phần bất
biến tìm được ở đây là dữ liệu đầu vào để giải trình hợp lý ở 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 47
Chương 6. GIẢI PHƯƠNG TRÌNH HỢP LÝ
Chương này đưa ra phương pháp giải phương trình hợp lý dựa vào tính bất biến
của cây sinh loài và mẫu dữ liệu quan sát.
6.1. Quỹ tích hợp lý trên một đa tạp
Chúng ta mô tả một mô hình thống kê là một tập con của:
1
0 1 0 1 0 1{( , ,..., ) : , ,..., 0 ... 1}
n
n n np p p p p p p p p
+Δ = ∈ > + + + =R vaø n
giả sử rằng, mô hình được mô tả như là một tập nghiệm chứa trong bởi một hệ các
phương trình các đa thức thuần nhất với các biến chưa biết . Các đa thức
được biết như là thành phần bất biến ở chương 5. Gọi V là tập của tất cả nghiệm phức
được cho bởi hệ các phương trình đa thức thuần nhất. Vấn đề cực đại hợp lý là tìm
những điểm ở mô hình
nΔ
0 1, ,..., np p p
0 1( , ,..., )np p p p=
0 nV V> = Δ∩
mà giải thích hợp lý nhất cho bởi véc tơ dữ liệu . Nghĩa là giải
quyết vấn đề tối ưu với ràng buộc sau:
1
0 1( , ,..., )
n
nu u u
+∈`
Cực đại hàm hợp lý 0 10 1 ... n
u u
n
uL p p p= hay bài toán log tương đương
với giả thuyết là 0 0log .... lognl u p u p= + + n 0p V>∈ .
Tiếp cận của chúng ta là tìm tất cả các điểm tới hạn của hàm hợp lý cực đại L
và sau đó chọn những nghiệm thực dương, những điểm đó là cực trị địa phương.
Trong quá trình giải tìm cực đại hàm trên, chúng ta sẽ tìm tất cả điểm tới hạn trên đa
tạp phức V. Cho sin gV ký hiệu những điểm kỳ dị của đa tạp V và tập .
Đặt P là iđêan trên vành đa thức được sinh bởi các đa thức được xác
định bởi V, hay:
sin: \reg gV V V=
0 1[ , ,..., ]np p pR
0 1[ ] [ , ,..., ]/nV p p p P=R R
Định nghĩa: Cho U là một tập mở 0 1\ ( . ... ( ))reg n iV p p p p℘ ∑ của V . Quỹ tích
hợp lý Zu là tập các điểm p U∈ mà 0dL = . Iđêan hợp lý [ ]uI V⊂ R là iđêan của tập
đóng của Zu trong V.
6.2. Ma trận Jacobi của các đa thức bất biến
6.2.1. Gradient- Vector vận tốc
Cho khả vi. Khi đó gradient của f tại x, được ký hiệu : nf →\ \ ( )gradf x và
định nghĩa là vector:
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 48
1
( ) ( ( ),..., ( ))grad
n
f ff x x
x x
x∂ ∂= ∂ ∂
Với , tập c∈\ 1{ : ( ) }nc ( )M x f x c f −= ∈ = =\ c gọi là mặt mức. Về mặt
hình học vector grad f(x) vuông góc với mặt mức của cM tại x.
Vậy phương trình mặt phẳng tiếp xúc với cM tại 1( ,..., )na a a= là
1 1
1
( )( ) ... ( )( ) 0n n
n
f fa x a a x a
x x
∂ ∂− + + − =∂ ∂
6.2.2. Ma trận Jacobi của các đa thức bất biến
Đặt là tập các đa thức thuần nhất được sinh ra bởi iđêan P.
Chúng ta mô tả ma trận Jacobi như sau:
0 1{ , ,..., }rg g g
1 0 1 1 1
2 0 2 1 2
0 1
1 1 1
n
n
r r r
g p g p g p
J g p g p g p
g p g p g p
⎡ ⎤⎢ ⎥∂ ∂ ∂ ∂ ∂ ∂⎢ ⎥⎢ ⎥= ∂ ∂ ∂ ∂ ∂ ∂⎢ ⎥⎢ ⎥⎢ ⎥∂ ∂ ∂ ∂ ∂ ∂⎣ ⎦
"
"
"
# # % #
" n
và ma trận chuyển vị của J là : TJ
1 0 2 0 0
1 1 2 1 1
1 2
1
1
1
r
rT
n n r
g p g p g p
g p g p g p
J
g p g p g p
∂ ∂ ∂ ∂ ∂ ∂⎡ ⎤⎢ ⎥∂ ∂ ∂ ∂ ∂ ∂⎢ ⎥= ⎢ ⎥⎢ ⎥∂ ∂ ∂ ∂ ∂ ∂⎣ ⎦
"
"
# # # % #
" n
Chúng ta nhân J bởi ma trận đường chéo các phần tử trên đường chéo là các
biến như sau:
0 1
1 1
0 1
0 1
2 2
0 10 1
0 1
0 1
0 0
. ( , ,..., )
n
n
n
nn
n
r r
n
n
p p p
g gp p p 1
2
r
g
p p p
g g gp p pJ J diag p p p
p p p
g g gp p p
p p p
⎡ ⎤⎢ ⎥∂ ∂ ∂⎢ ⎥∂ ∂ ∂⎢ ⎥⎢ ⎥∂ ∂ ∂⎢ ⎥= = ∂ ∂ ∂⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥∂ ∂ ∂⎢ ⎥∂ ∂ ∂⎣ ⎦
"
"
"
# # % #
"
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 49
6.2.3. Không gian tiếp xúc
Vì V là một đa tạp khả vi và p V∈ , tập vector tiếp xúc với M tại p được gọi là
không gian tiếp xúc với V tại p và ký hiệu thì: pT V
{ : ( ), 1,...,gradnp iT V v v g p i m= ∈ ⊥ =P }
Viết một các khác cho bởi hệ phương trình pT V : . 0
nv J v∈ =P , hay là
Ker của ma trận trên .
pT V
J 0 1[ ] [ , ,..., ]/nV p p p=\ R P
6.3. Bài toán cực trị điều kiện
Ta có bài toán tìm cực trị của hàm hợp lý chính là tìm cực trị hàm .
Nói cách khác là tìm cực trị của hàm l với điều kiện ràng buộc .
Ta có:
:l V →\
0 1 ... 0rg g g= = = =
0
0
( ) ( ,..., )grad n
n
u ul p
p p
=
Điều kiện cần: Nếu đặt cực trị với ràng buộc l 0 1 ... 0rg g g= = = = tại p, thì
( )grad pf p T V⊥ suy ra
0
0
0
... 0n n
n
u uv v
p p
+ + = với pv T V∈
vì nên tập các p thuộc tập: 0, 0,...,ip i≠ = n
0
0
n
i i
i
uφ
=
=∑ với các vector 0( ,..., )nφ φ
chạy trên tập Ker của ma trận trên J 0 1[ ] [ , ,..., ]/nV p p p P=\ R
Phương pháp nhân tử hoá Lagrange: Từ trên, để tìm điểm nghi ngờ cực trị của l với
điều kiện , ta lập hàm Lagrange 0 1 ... 0rg g g= = = =
1
0 0 0( , ) ( ) ... , ,( ,..., )
r
r r rp l p g g p Vλ λ λ λ λ += − − − ∈ ∈^L
Nếu p là điểm cực trị có điều kiện thì tồn tại sao cho ( ,10( ,..., )
r
rλ λ +∈^ )p λ là
nghiệm của hệ:
0
( , ) 0
( ) 0
( ) 0r
p
p
g p
g p
λ∂⎧ =⎪ ∂⎪⎪ =⎨⎪⎪ =⎪⎩
#
L
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 50
Điều kiện đủ: Đặt *( *, )pH p λL là Hessian của hàm Lagrange theo biến p. Khi
đó: Nếu
L
*( *, )pH p λL xác định âm, thì đạt cực đại tại . l *p
( Xem [Tạ Lê Lợi, 2002])
6.4. Bậc của hợp lý cực đại
Bậc hợp lý cực đại của mô hình thống kê đại số là số điểm tới hạn phức của
logarit phương trình hợp lý đối với vector dữ liệu l 1nu +∈` và nó hữu hạn và bị chặn
trên theo định lý sau:
Định lý: Cho 0,..., rP g g= 〈 〉 là một iđêan trong với là đa thức thuần
nhất có bậc là với . Thì bậc hợp lý cực đại của P là hữu hạn và bị chặn
trên bởi
0[ ,..., ]np p_ ig
id 0,...,i = r
d0 1
0 1
0
0 0
...
0,...., 0
... r
r
r
i i i
r
i i i n r
i i
D d d
+ + + ≤ −
> >
= ∑
Chứng minh: Xem [Hoşten et al., 2005].
6.5. Các thuật toán
Từ các phân tích trên ta xây dựng một số thuật toán sau giải phương trình hợp
lý.
Thuật toán 1: (Tính toán phương trình hợp lý)
Input: Một iđêan thuần nhất và một vector 0[ ,..., ]nP p p⊂ \ 1nu +∈` .
Output: Iđêan hợp lý uI của mô hình ( )V Pυ= cho dữ liệu u.
Bước 1: Tính . Đặt Q là iđêan những quỹ tích điểm kỳ dị của V. ( 1) dim(c n V= + − )
Bước 2: Tính Ker M của ma trận trên J 0 1[ ] [ , ,..., ]/nV p p p P=\ R .
Bước 3: Đặt 'uI là iđêan trong sinh bởi đa thức [ ]V\
0
0
n
i i
i
uφ
=
=∑ , với vector
0( ,..., )nφ φ chạy trên tập M.
Bước 4: Iđêan uI bằng
'
uI loại bỏ đi những điểm kỳ dị.
(Chứng minh tính đúng thuật toán trên: Bước 1, tìm quỹ tích những điểm kỳ dị được
chứng minh ở [J.S. Milne, 2005]. Các bước còn lại có được bởi sự phân tích ở phần
đầu chươ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 51
Thuật toán 2: (Tính toán cực đại địa phương của hàm hợp lý)
Input: Iđêan hợp lý uI của mô hình V và dữ liệu u.
Output: Danh sách của tất cả các cực đại địa phương cho phương trình hợp lý.
Bước 1: Nếu dim đối với dữ liệu u, tính ra tập ( ) 0uI = uZ .
Với mỗi một nghiệm dương thực hiện các bước sau: * 0up Z V>∈ ∩
Bước 2: Giải phương trình tuyến tính *( ).TJ p uλ = để thu các nhân tử Lagrange *iλ .
Bước 3: Tính toán Hessian *( *, )pH p λL của hàm ( , )p λL .
Bước 4: Nếu *( *, )pH p λL ở bước 3 là xác định âm thì xuất với giá trị hàm
tương ứng.
*p
*(l p )
6.6. Áp dụng giải phương trình hợp lý
Xét ví dụ cây sinh loài ở chương 5, thành phần bất biến chúng ta tìm được là:
2 2 2 3 2 2
1 0 2 0 1 0 2 1 1 2 1 2 2
8 16 4 80 8 4 4- - -
3 27 9 729 81 27 27
g p p p p p p p p p p p p= + + + 3 0=
và bất biến tầm thường là:
0 0 1 2 1 0g p p p= + + − =
hay :
2 2 2 3 2 2
0 1 2 0 2 0 1 0 2 1 1 2 1 2 2
8 16 4 80 8 4 4: 1, - - -
3 27 9 729 81 27 27u
3I p p p p p p p p p p p p p p p= 〈 + + − + + + 〉
Giả sử bộ dữ liệu quan sát đếm được cho từng xác suất tương ứng là
.
0 1 2, ,p p p
(700,45,49)u =
Kết quả sau khi chạy thuật toán 1: Xem phụ lục 2.
Kết quả sau khi chạy thuật toán 2 tìm được số nghiệm của phương trình hợp lý
là 5 nghiệm như sau:
(0.2474737682, 2.5624447362, -1.8099185045)
(0.8289181084, 0.1641242871, 0.0069576045)
(0.1469015638, 0.7993903746, 0.0537080615)
(-0.0296031392-i*0.2337477875,0.0376058385-i*0.0038708209,
0.9919973007+i*0.2376186084)
(-0.0296031392+i*0.2337477875,0.0376058385+i*0.0038708209,
0.9919973007-i*0.2376186084)
Trong 5 nghiệm trên có 2 nghiệm thực dương là:
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 52
(0.8289181084, 0.1641242871, 0.006957604546)
(0.1469015638, 0.7993903746, 0.05370806152)
Hàm Lagrange tương ứng:
0 1 2 0 0( , ) 700ln( ) 45ln( ) 49ln( )p p p p g 1 1gλ λ λ= + + + +L
Với nghiệm thứ nhất nhân tử Lagrange tìm được tương ứng là:
*
0 844.4742431640λ = , *1 =18.6903667449λ
Ma trận Hessian tương ứng
* *
0 1 2
1( , )
-7.000e+02 -3.441e-03 7.832e-02
-3.441e-03 -4.502e+01 -6.100e-04
7.832e-02 -6.100e-04 -4.899e+01
pH p p p p
λ
⎡ ⎤⎢ ⎥= ⎢ ⎥⎢ ⎥⎣ ⎦
L
Ma trận Hessian trên cho xác định âm, cho nên hàm hợp lý trên đạt giá trị cực
đại cục bộ tại điểm trên và giá trị hàm hợp lý tương ứng là:
-456.0927286l =
Với nghiệm thứ hai nhân tử Lagrange tìm được tương ứng là:
*
0 4765.0976562500λ = , *1 =64.9804763793λ
Ma trận Hessian tương ứng
* *
0 1 2
1( , )
-6.999e+02 -3.883e-01 3.407e-01
-3.883e-01 -4.486e+01 -7.124e-02
3.407e-01 -7.124e-02 -4.904e+01
pH p p p p
λ
⎡ ⎤⎢ ⎥= ⎢ ⎥⎢ ⎥⎣ ⎦
L
Ma trận Hessian trên cũng cho xác định âm, cho nên hàm hợp lý trên đạt giá trị
cực đại cục bộ tại điểm trên và giá trị hàm hợp lý tương ứng là:
-1495.955966l =
Với 2 kết quả trên ta có thể kết luận, hàm hợp lý đạt giá trị cực đại toàn cục tại
(0.8289181084, 0.1641242871, 0.0069576045) p = và giá trị hàm hợp lý tương
ứng là . -456.0927286l =
Kết quả chương trình với số cây khác xem phụ lục 4.
Vì chúng ta xử lý bài toán trên những xác suất nên độ chính xác là cần thiết, vì
thế độ dài phần lẻ được lấy ở chương trình là 10 con số. Và chương 7 sẽ cho chúng ta
một cái nhìn tổng quan hơn về chương trình thực hiện giải phương trình hợp lý - áp
dụng 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 53
Chương 7. CHƯƠNG TRÌNH THỰC HIỆN
7.1. Sơ đồ khối chương trình
Hình 9: Sơ đồ khối chương trình tìm cấu trúc cây sinh loài
Input: Các chuỗi DNA
N := Số cây sinh loài được
tạo ra từ các chuỗi trên
i: = 1
Tính các thành phần bất biến
của cây i
Tính vector u từ các chuỗi
Begin
Giải phương trình hợp lý
i<N
Output: Cây sinh loài có hàm hợp
lý lớn nhất
Tính các lớp tương đương trên
cây i
End
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 54
7.2. Sơ lược về chương trình
Bài toán trên được chia thành nhiều bài toán nhỏ tách rời nhau. Vì thế để chọn
một ngôn ngữ lập trình phù hợp của từng công đoạn cũng cần được tính đến. Trong
luận văn này, hai ngôn ngữ được chọn để viết và kiểm nghiệm đó là ngôn ngữ C++ và
Singular. Gói phần mềm viết bằng C++ giải quyết phần đầu bài toán là xác định các
xác suất xảy ra trên cây và đếm số lượng các vector dữ liệu u.
Phần quan trọng nhất của chương trình là tìm các thành phần bất biến và giải
phương trình hợp lý được viết trên ngôn ngữ Singular (Xem phụ lục 3). Singular là một
hệ thống các phương pháp tính toán đa thức, đại số giao hoán và hình học đại số. Nó là
một gói phần mềm rất hữu hiệu cho đại số thống kê.
7.3. Kết quả chương trình
Xét cây sinh loài có gốc với 3 taxa, dữ liệu tương ứng với mỗi lá ở đây là những
gen có tên là HIVenvSweden (dữ liệu lấy ở [Yang, Z. 1997]) có chiều dài 273 ký tự
như sau:
U68496:=GTAGTAATTAGATCTGAAAACTTCACGAACAATGCTAAAACCAT
AATAGTACAGCTAAATAAATCTGTAGAAATTAATTGTGTAAGACCCGGCA
ACAATACAAGAAGAAGTATACATATAGGACCAGGGAGAGCATATTATACA
GGAGAAGTAATAGGAGATATAAGACAAGCACATTGTAACCTTAGTAGAAC
AGACTGGAATAAAACTTTAAAACAGGTAGCTGAAAAATTAAGAGAACAAT
TTAATACAACAATAGTCTTTAATCAATCC
U68497:=GTAGTAATTAGATCTGAAAACTTCTCGAACAATGCTAAAACCATA
ATAGTACAGCTAAATAAATCTGTAGAAATTAATTGTACAAGACCCAACAA
CAATACAAGAAGAAGTATACATTTTGGACCAGGGAAAGCATTTTATGCAG
GAGAAATAATAGGAGATATAAGACAAGCATATTGTACCCTTAATGGAACA
GAATGGAATAACACTTTAAAACAGGTAGCTGAAAAATTAAGAGAACAATT
TATTAAAACAATAGTTTTTAATCAATCC
U68498:=ATAGTAATTAGATCTGAAAACTTCTCGAACAATGCTAAAACCATA
ATAGTACAGCTAAATAAATCTGTAGAAATTAATTGTACAAGACCCAACAA
CAATACAAGAAGAAGTATACATTTCGGACCAGGGAAAGCATTTTATGCAG
GAGAAATAATAGGAGATATAAGACAAGCATATTGTACTCTTAATGGAGCA
GAATGGAATAACACTGTAAAACAGGTAGCTGCAAAATTAAGAGAACAATT
TAATAAAACAATAATCTTTAATCAATCC
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 55
Với 3 gen này mục đích chúng ta là tìm mối liên hệ tổ tiên giữa chúng. Có 2
dạng cây 3 taxa có gốc như hình 10.
Cây 3 taxa hình móng Cây3 taxa hình lược
Hình 10: Hai hình dạng cây 3 taxa có gốc
Ta chọn mô hình Jukes – Cantor cho bài toán này.
Với Cây hình móng: Chỉ có một trường hợp duy nhất, dữ liệu quan sát khi so
dóng trên từng cột như Bảng 2.
Bảng 2: Các mẫu và số lượng từng mẫu trên 3 chuỗi gen HIVenvSweden với cây hình
móng (U68496, U68497, U68498)
Các
mẫu
A
A
A
A
A
T
A
G
A
A
G
G
A
C
A
A
C
C
A
C
T
A
T
T
G
A
A
G
A
G
G
G
G
C
A
A
C
C
C
C
C
T
C
T
C
C
T
T
T
G
T
T
C
C
T
T
T
Số
lượng 113 1 1 2 1 2 1 3 6 2 39 2 34 1 1 1 1 1 61
Cây này có 3 lớp tương đương, số lượng từng lớp trên dữ liệu là u=(247, 25, 1),
khi chạy chương trình cho kết quả:
p= (0.903285626594, 0.094557322708, 0.00215705069738)
và giá trị hàm hợp lý tương ứng ( ) -90.22670674l p =
Với Cây hình lược: Có 4 lớp xác suất tương đương và 3 trường hợp xảy ra:
Trường hợp 1: Các nút cây tương ứng (U68496,(U68496, U68496)) và số
lượng các mẫu dữ liệu quan sát khi so dóng trên từng cột đối với trường hợp này ở
bảng 3.
Bảng 3: Các mẫu và số lượng từng mẫu trên 3 chuỗi gen HIVenvSweden với cây hình
lược với trường hợp ((U68496,(U68497, U68498))
Các
mẫu
A
A
A
A
A
T
A
G
A
A
G
G
A
C
A
A
C
C
A
C
T
A
T
T
G
A
A
G
A
G
G
G
G
C
A
A
C
C
C
C
C
T
C
T
C
C
T
T
T
G
T
T
C
C
T
T
T
Số
lượng 113 1 1 2 1 2 1 3 6 2 39 2 34 1 1 1 1 1 61
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 56
Vector dữ liệu tương ứng là: u = (247, 8, 17, 1) và kết quả khi chạy chương
trình:
p =(0.9027837320165, 0.03154973251390, 0.0640321523858, 0.001634383083735)
và giá trị hàm ( ) -106.0495469l p =
Trường hợp 2: Các nút cây tương ứng (U68497,(U68496, U68498)) và số
lượng các mẫu dữ liệu quan sát khi so dóng trên từng cột đối với trường hợp này giống
như bảng 2, vì thế tương tự như trường hợp 1.
Trường hợp 3: Các nút cây tương ứng (U68498,(U68496, U68497)) và số
lượng các mẫu dữ liệu quan sát khi so dóng trên từng cột đối với trường hợp này trên
bảng 4.
Bảng 4: Các mẫu và số lượng từng mẫu trên 3 chuỗi gen HIVenvSweden với cây hình
lược với trường hợp ((U68498,(U68496, U68497))
Các
mẫu
A
A
A
A
A
G
A
A
C
A
G
G
A
T
A
G
A
A
G
G
A
G
G
G
G
G
T
C
A
A
C
C
A
C
C
C
C
C
T
C
T
A
C
T
C
T
C
C
T
T
A
T
T
C
T
T
T
Số
lượng 113 6 2 2 1 1 2 39 1 1 2 34 1 1 1 1 3 1 61
Vector dữ liệu tương ứng u = (247, 19, 6, 1), kết quả khi chạy chương trình:
p =(0.9032484236332, 0.07179349581611, 0.0228400591792, 0.0021180213714)
và giá trị hàm ( ) -104.0121158l p =
Với kết quả trên thì ta thấy với trường hợp cây sinh loài thứ 2 thì các dữ liệu
trên đựợc gán cho cây với thứ tự ((U68498,(U68496, U68497)) là “hợp lý” nhất. Xét
về cây thì cây sinh loài giải thích tốt nhất cho dữ liệu trên là cây sinh loài hình lược, vì
nó có hàm hợp lý cho kết quả lớn nhấ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 57
Chương 8. TỔNG KẾT – ĐÁNH GIÁ
Chương này tổng kết lại những công việc đã làm được, sau đó nêu ra những
đóng góp và hướng phát triển của luận văn.
8.1. Tổng kết
Luận văn đã nghiên cứu đến một lãnh vực đang được rất nhiều người quan tâm
trong giai đoạn hiện nay. Luận văn cho ta một cái nhìn tổng thể về cây sinh loài và
phương pháp để giải phương trình hợp lý, cũng như cách mô hình hóa bài toán và chọn
phương pháp và giải quyết bài toán. Trong đó phương pháp được chọn là phương pháp
đại số nói chung và đại số máy tính nói riêng. Ngoài ra, luận văn cũng cung cấp một
gói phần mềm được sử dụng như là kiểm nghiệm những lý thuyết đưa ra.
Vì kiến thức cũng như thời gian có hạn nên luận văn có nhiều khiếm khuyết
như nội dung luận văn chưa được cô đọng, phần chương trình kiểm nghiệm các mô
đun chưa tích hợp được các phần lại với nhau. Hạn chế lớn nhất phải kể đến là chương
trình chỉ chạy được trên cây sinh loài 3 taxa, còn trường hợp 4 và 5 taxa chỉ chạy được
ở những mô hình đơn giản.
8.2. Những đóng góp của luận văn
- Đối với bản thân qua luận văn này đã thu thập kiến thức cơ bản về sinh học
nói chung và việc nghiên cứu tiến hoá hiện nay nói riêng. Cũng qua luận văn này đã
bổ sung cho tôi kiến thức về toán học, xác suất thống kê và nhất là đại số máy tính.
- Luận văn này như là một tài liệu tham khảo ban đầu cho những ai quan tâm
đến lãnh vực tính toán thống kê về sinh học phân tử, về cây tiến hóa. Ngoài ra cũng là
tài liệu bổ ích cho một cái nhìn tổng quan về phương pháp ước lượng hợp lý cực đại và
cách giải chúng bằng phương pháp đại số. Những lãnh vực được đề cập trong luận văn
này đang được quan tâm và phát triển mạnh trên thế giới trong những năm gần đây.
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 58
8.3. Hướng phát triển
Với những hạn chế đã nêu ở trên, hướng phát triển luận văn đưa ra chủ yếu là
khắc phục những nhược điểm là làm sao để giải bài toán trên những cây sinh loài lớn
hơn:
- Cải tiến các thuật toán để có thể chạy tốt hơn.
- Song song hóa để giải bài toán ước lượng hợp lý cực đại.
- Một hướng khác cũng nên quan tâm là nghiên cứu đặc tính của từng cây cụ thể
và tập trung giải quyết bài toán trên những cây riêng biệ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 59
TÀI LIỆU THAM KHẢO
[Hoàng Xuân Sính,
2003]
Hoàng Xuân Sinh. Đại số đại cương. Nhà xuất bản GD, 2003.
[Đào Hữu Hồ, 2002] Đào Hữu Hồ. Xác suất thống kê. Nhà xuất KHKT, 2002.
[Tạ Lê Lợi, 2002] Tạ Lê Lợi. Giải Tích – Giáo trình cử nhân toán. Đại học Đà Lạt, 2002.
[Benny Chor, 2005] Maximum likelihood Jukes-Cantor Triplets: Analytic Solutions. arXiv:
q-bio.PE/0505054 v1 27 May 2005
[B. Chor, M. Hendy, B.
Holland, and D. Penny,
2000]
Multiple maxima of likelihood in phylogenetic trees: An analytic
approach. MBE, 17(10):1529–1541, 2000. Earlier version appeared in
RECOMB 2000.
[Catanese et al., 2005]
F Catanese, S Ho¸sten, A Khetan, and B Sturmfels. The maximum
likelihood degree. American Journal of Mathematics, 2005. To appear.
Cited on p. 10, 106
[Evans and
Speed,1993]
Steven N.Evans and Speed. Invariants of some probability models used
in phylogenetic inference. The Annals of Statistics 1993, vol. 21, No.
1, 355-377.
[Hoşten et al., 2005]
S Hoşten, A Khetan, and B Sturmfels. Solving the likelihood
equations. Foundations of Computational Mathematics, 2005. To
appear. Cited on p. 108, 335, 343, 344
[J.S. Milne, 2005] J.S. Milne. Algebraic Geometry. Taiaroa Publishing Erehwon. Version
5.00. February 20,2005
[Lior Pachter and
Bernd Sturmfels, 2005]
Lior Pachter and Bernd Sturmfels. Algebraic Statistics for
Computational Biology. Cambridge University Press 2005.
[Luis D. Garcia, 2006] Luis David Garcia. Small Phylogennetic Trees.
[Melvyn B.
Nathanson,1999]
Melvyn B. Nathanson. Elementary Methods in Number Theory .
Springer, 1999.
[Nguyen V.M.
Man,2007]
Nguyen V.M. Man. Computational Algebraic Statistics(CAS).
[Sturmfels and
Sullivant, 2005]
Bernd Sturmfels and Seth Sullivant. Toric Ideals of Phylogenetic
Invariants.
[Yang, Z. 1997] Ziheng Yang. A program package for phylogenetic analysis by
maximum likelihood. Computer Applications in BioSciences 13:555-
556
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 60
Phụ lục 1. Tập các xác suất trình bày ở
chương 5
pAAA=1/4*c0*c0*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^3+3/4*c1^3,
pAAG=1/4*c0*c1*c0+1/4*c1*c0*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pAAC=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pAAT=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c1*c0*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pAGA=1/4*c1*c0*c0+1/4*c0*c1*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pAGG=1/4*c1*c1*c0+1/4*c0*c0*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pAGC=1/4*c1*c1*c0+1/4*c0*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pAGT=1/4*c1*c1*c0+1/4*c0*c1*c1+1/4*c1*c1*c1+1/4*c1*c0*c1=3/4*c0*c1^2+1/4*c1^3,
pACA=1/4*c1*c0*c0+1/4*c1*c1*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pACG=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pACC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c0*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pACT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c0*c1*c1+1/4*c1*c0*c1=3/4*c0*c1^2+1/4*c1^3,
pATA=1/4*c1*c0*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c0*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pATG=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c1*c1*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pATC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pATT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c0*c0*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGAA=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGAG=1/4*c0*c1*c0+1/4*c1*c0*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGAC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pGAT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pGGA=1/4*c1*c0*c0+1/4*c0*c1*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGGG=1/4*c0*c0*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^3+3/4*c1^3,
pGGC=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGGT=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGCA=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pGCG=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGCC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c0*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGCT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c0*c1*c1+1/4*c1*c0*c1=3/4*c0*c1^2+1/4*c1^3,
pGTA=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pGTG=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pGTC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pGTT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c0*c0*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCAA=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCAG=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pCAC=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCAT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
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 61
pCGA=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pCGG=1/4*c1*c1*c0+1/4*c0*c0*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCGC=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCGT=1/4*c1*c1*c0+1/4*c0*c1*c1+1/4*c1*c1*c1+1/4*c1*c0*c1=3/4*c0*c1^2+1/4*c1^3,
pCCA=1/4*c1*c0*c0+1/4*c1*c1*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCCG=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCCC=1/4*c0*c0*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^3+3/4*c1^3,
pCCT=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCTA=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pCTG=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c1*c1*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pCTC=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pCTT=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c0*c0*c1=1/4c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTAA=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTAG=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pTAC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c0*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pTAT=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c1*c0*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTGA=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pTGG=1/4*c1*c1*c0+1/4*c0*c0*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTGC=1/4*c1*c1*c0+1/4*c0*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pTGT=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTCA=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pTCG=1/4*c1*c1*c0+1/4*c1*c0*c1+1/4*c0*c1*c1+1/4*c1*c1*c1=3/4*c0*c1^2+1/4*c1^3,
pTCC=1/4*c1*c1*c0+1/4*c1*c1*c1+1/4*c0*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTCT=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTTA=1/4*c1*c0*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c0*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTTG=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTTC=1/4*c0*c1*c0+1/4*c1*c1*c1+1/4*c1*c0*c1+1/4*c1*c1*c1=1/4*c0^2*c1+1/4*c0*c1^2+1/2*c1^3,
pTTT=1/4*c0*c0*c0+1/4*c1*c1*c1+1/4*c1*c1*c1+1/4*c1*c1*c1=1/4*c0^3+3/4*c1^3;
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 62
Phụ lục 2. Tập các dữ liệu kết quả thực hiện
trình bày ở chương 6
c = (n+1) - dim(V) = 2
Q:= <-192*p1^2 - 576*p1*p2 - 432*p2^2+756*p1+810*p2 - 486,
- 576*p1^2 - 1728*p1*p2 - 1296*p2^2+972*p1+1782*p2 - 486,
384*p1^2+1152*p1*p2+864*p2^2 - 216*p1 - 972*p2>
M = Ker( ) := <775168*p1^2*p2^2*gen(3) +11627520*p1^2*p2^2*gen(2) +
204649114632192*p1^2*p2^2*gen(1) + 2325504*p1*p2^3*gen(3) + 66276864*p1*p2^3*gen(2) +
767434179870720*p1*p2^3*gen(1) + 1744128*p2^4*gen(3)+52323840*p2^4*gen(2) +
920921015844864*p2^4*gen(1) + 1975872*p1^2*p2*gen(3) + 14539584*p1^2*p2*gen(2) +
3700197350079168*p1^2*p2*gen(1) - 1539648*p1*p2^2*gen(3) - 50959872*p1*p2^2*gen(2) +
1839387248964864*p1*p2^2*gen(1) - 3920832*p2^3*gen(3) - 101259504*p2^3*gen(2) -
1822005910233648*p2^3*gen(1) - 2102400*p1^2*gen(3) - 60298560*p1^2*gen(2) -
2339638635463488*p1^2*gen(1) - 2890566*p1*p2*gen(3) - 69060222*p1*p2*gen(2) -
5046845491193202*p1*p2*gen(1) + 2629503*p2^2*gen(3) + 45131337*p2^2*gen(2) +
876411667623555*p2^2*gen(1) + 1773900*p1*gen(3) + 50876910*p1*gen(2) +
1974070098672318*p1*gen(1) - 364014*p2*gen(3) + 7490718*p2*gen(2) + 87067895564682*p2*gen(1),
775168*p1^2*p2^2*gen(3) + 10077184*p1^2*p2^2*gen(2) + 2325504*p1*p2^3*gen(3) +
70927872*p1*p2^3*gen(2) + 51162278658048*p1*p2^3*gen(1) + 1744128*p2^4*gen(3) +
62788608*p2^4*gen(2) + 306973671948288*p2^4*gen(1) + 3518144*p1^2*p2*gen(3) +
37749568*p1^2*p2*gen(2) + 5725315696773056*p1^2*p2*gen(1) + 1246080*p1*p2^2*gen(3) +
3243264*p1*p2^2*gen(2) + 6256978794324480*p1*p2^2*gen(1) - 2340288*p2^3*gen(3) -
80025840*p2^3*gen(2) + 1456573961040912*p2^3*gen(1) - 3132672*p1^2*gen(3) - 85125056*p1^2*gen(2) -
3356546566796480*p1^2*gen(1) - 8246130*p1*p2*gen(3) - 151489470*p1*p2*gen(2) -
12422421565374378*p1*p2*gen(1) - 1598535*p2^2*gen(3) - 27016119*p2^2*gen(2) -
5667479873380377*p2^2*gen(1) + 4745592*p1*gen(3) + 132122826*p1*gen(2) +
5171724801198018*p1*gen(1) + 4078296*p2*gen(3) + 93163662*p2*gen(2) + 6095180048475222*p2*gen(1)
- 1773900*gen(3) - 50876910*gen(2) - 1974070098672318*gen(1)
J
, 775168*p1^2*p2^2*gen(3) + 2325504*p1*p2^3*gen(3) + 45347328*p1*p2^3*gen(2) +
1744128*p2^4*gen(3) + 47091456*p2^4*gen(2) + 230230253961216*p2^4*gen(1) +
5355200*p1^2*p2*gen(3) + 54461952*p1^2*p2*gen(2) + 8363661940680768*p1^2*p2*gen(1) +
4964640*p1*p2^2*gen(3) + 29969856*p1*p2^2*gen(2) + 13016976252241536*p1*p2^2*gen(1) -
459792*p2^3*gen(3) - 64457856*p2^3*gen(2) + 5516397190398384*p2^3*gen(1) - 4660224*p1^2*gen(3) -
123711552*p1^2*gen(2) - 4949039207382912*p1^2*gen(1) - 15788844*p1*p2*gen(3) -
289782558*p1*p2*gen(2) - 22584712180858542*p1*p2*gen(1) - 7662492*p2^2*gen(3) -
104848263*p2^2*gen(2) - 15607396038950235*p2^2*gen(1) + 8045784*p1*gen(3) + 222790392*p1*gen(2) +
8734948040514888*p1*gen(1) + 10033956*p2*gen(3) + 223300476*p2*gen(2) +
13978410003807732*p2*gen(1) - 3547800*gen(3) - 101753820*gen(2) - 3948140197344636*gen(1)
'
uI :=<12013588416*p2^3 + 148826033904*p1^2 + 118484015288*p1*p2 - 9505114572*p2^2 -
149625916917*p1 - 20725910325*p2 + 20557713033, 3325964017042946421278092794380392608*p1*p2^2
+ 4983618577401699676273006787935157328*p2^3 + 22774007877394191975227302283935849888*p1^2 +
14264618759069671514915176286636353554*p1*p2 - 7849236572356881083751999247248463887*p2^2 -
22881274699550634255631609142426367189*p1 + 293064477807175251834460143692430423*p2 +
3123938321374671114025459142247553701, 296521570634524577657257064180480*p1^2*p2 +
678645469410592180933991906124096*p1*p2^2 + 399557794965576369841442093252256*p2^3 +
163472584389053652220661632656096*p1^2 - 626014046035479156390619894451652*p1*p2 -
805886945003862357039771874801974*p2^2 - 161563250873918603863505188608723*p1 +
412431155271902586749997779755821*p2 + 19936260526014638257361883471927>
uI :=<57425888*p1^2 + 60331008*p1*p2 + 5850720*p2^2 - 58311573*p1 - 13781205*p2 + 8050185,
43118186428457088*p2^3 - 135923428146969968*p1*p2 - 88536193620103176*p2^2 +
5367437095038981*p1 + 53799892346655045*p2 – 1095761223392121, 28745457618971392*p1*p2^2 +
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 63
52276369682092096*p1*p2 + 549084854058576*p2^2 - 3253058199626023*p1 - 3973831866790983*p2 +
501604238561283, p0+p1+p2-1,128*p1^3 + 576*p1^2*p2 + 864*p1*p2^2 + 432*p2^3 - 108*p1^2 -
972*p1*p2 - 891*p2^2 + 486*p2 + 19936260526014638257361883471927>
Trong đó: gen(n) =[0,0,….,1], ví dụ gen(3)=[0,0,1].
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 64
Phụ lục 3. Trích một số SourceCodes chương
trình viết trên Singular
//Chương trình giải quyết cây 3 taxa hình móng.
LIB "control.lib";
LIB "solve.lib";
killall();
int i,j;
//Đây là các bất biến q của phép biến đổi Fourier.
ring rQ = 0,(q1,q2,q3),dp;
ideal Invariants =
q2^3-q1*q3^2;
// Đây là các lớp xác suất tương đương p.
ring rP = 0,(p1,p2,p3),dp;
//Ma trận chuyển đổi Fourier.
matrix qtop[3][3] =
1,1,1,
1,1/9,-1/3,
1,-1/3,1/3;
ideal Fourier = qtop*transpose(maxideal(1));
// Danh sách các đa thức bất biến
map F = rQ, Fourier;
ideal PInvariants = F(Invariants);
ideal V=p1+p2+p3-1,PInvariants;
matrix J=jacob(V);
J;
"dim(V)=";dim(V);
ideal P=groebner(V);
P;
"dim(P)=";dim(P);
// codim của V
int c=nvars(rP)-dim(P);
"c=";c;
matrix J=jacob(P);
ideal Q=minor(J,c);
"Q=";print(Q);
matrix Jn=J*transpose(diag(maxideal(1)));
"Jnt=";Jn;
qring q=std(P);
matrix Jn=fetch(rP,Jn);
"Jns=";Jn;
//M là nhân của ma trận chuyển vị Jn
module M=transpose(rightKernel(Jn));
"M=";M;
//Vector dữ liệu u
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 65
matrix u[3][1]=247,25,1;
ideal Ipu=groebner(flatten(M*u));
"Ipu=";print(Ipu);
dim(Ipu);
Ipu;
ring rr = 0, (p1,p2,p3), dp;
ideal Ipu=fetch(q,Ipu);
ideal Q=fetch(rP,Q);
//I chứa các nghiệm của phương trình hợp lý sau khi loại các điểm kỳ dị
list I=sat(Ipu,p1*p2*p3*(p1+p2+p3));
I;dim(I[1]);
for(i=2;i<=size(Q);i++)
{
list A=sat(I[1],p1*p2*p3*(p1+p2+p3)*Q[i]);
A;
I=A;
"dim(I)=";dim(I[1]);vdim(I[1]);
}
ideal X=I[1],fetch(rP,P);
dim(X);
X;
def AC=solve(X,30,0,"nodisplay");
setring AC;
//SOL chứa danh sách các nghiệm hợp lý
SOL;
matrix Jn=transpose(fetch(rP,Jn));
Jn;
matrix u=fetch(q,u);
u;
//Hàm trả về ma trận J sau khi thay thế các nghiệm vào
proc Replace(matrix M,list SOL)
{
matrix Y=M;
int ii,jj,kk;
for(ii=1;ii<=nrows(M);ii++)
{
for(jj=1;jj<=ncols(M);jj++)
{
for(kk=1;kk<=nvars(rP);kk++)
{
Y[ii,jj]=subst(Y[ii,jj],var(kk),SOL[kk]);
}
}
}
return (Y);
}
matrix Y;
//Đoạn chương trình chọn những nghiệm thực dương phương trình hợp lý
for(i=1;i<=size(SOL);i++)
{
for(j=1;j<=size(SOL[i]);j++)
{
if(SOL[i][j]=1)
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 66
{
break;
}
}
if(j>size(SOL[i]))
{
SOL[i];
matrix Y=Replace(Jn,SOL[i]);
Y;
}
}
//Đoạn chương trình tính các nhân tử Lagrange
LIB "presolve.lib";
LIB "solve.lib";
ring r = real, (x,y), lp;
matrix Y[4][2];
//Chèn file ma trận Y tính được ở trên
matrix K=Y*transpose(diag(maxideal(1)));
K;
matrix A[2][1]=1,1;
matrix I=K*A;
I;
ideal M=flatten(I);
M;
ideal j= 113,46,114;
j = M-j;
ideal z=solvelinearpart(j,0);
z;
def AC=solve(z,30,0,"nodisplay");
setring AC;
SOL;
//Xác định ma trận Hessian là xác đinh âm hay không?
LIB "linalg.lib";
ring r=real,(p1,p2,p3),lp;
//Chèn file poly f bằng nhân tử Lagrange nhân với các gi
matrix J=jacob(jacob(f));
print(J);
//Ma trận I la một nghiệm của p
//matrix I[3][1]=0.09439863069,0.6214192685,0.2841821009;
I;
int u,v,k;
for(u=1;u<=nrows(J);u++)
{
for(v=1;v<=ncols(J);v++)
{
for(k=1;k<=nvars(rP);k++)
{
J[u,v]=subst(J[u,v],var(k),I[k,1]);
}
}
}
print(J);
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 67
" ";
for(u=1;u<=nrows(J);u++)
{
for(v=1;v<=ncols(J);v++)
{
for(k=1;k<=nvars(rP);k++)
{
J[u,v]=J[u,v]*I[k,1];
}
}
}
print(J);
" ";
//vector u
//matrix l[3][1]=-700,-45,-49;
matrix H[3][3]=diag(l)-J;
print(H);
//Xác định âm ma trận Hessian H, nếu âm trả về 0
pos_def(H);
Một số hàm quan trọng sử dụng trong chương trình
minor(matrix_expression,int_expression):returns the set of all minors (=subdeterminants) of
the given size of a matrix.
rightKernel(M): computes the right kernel of matrix M (a module of all elements v such that
Mv=0).
jacob( ideal_expression ):computes the Jacobi ideal, resp. Jacobi matrix, generated by all
partial derivatives of the input.
sat(id,j); id=ideal/module, j=ideal: return list of an ideal/module [1] and an integer [2]: [1] =
saturation of id with respect to j (= union_(k=1...) of id:j^k) [2] = saturation exponent (= min(
k | id:j^k = id:j^(k+1) ))
solve(G [, m, n [, l]] [,"oldring"] [,"nodisplay"] ); G = ideal, m, n, l = integers (control
parameters of the method), outR ring m: precision of output in digits ( 4 <= m) and of the
generated ring of complex numbers; n: control of multiplicity n = 0, return all different roots n
!= 0, find all roots (with multiplicity) l: precision of internal computation in decimal digits ( l
>=8 ) only if the basering is not complex or complex with smaller precision.
solvelinearpart(id [,n] ); id=ideal/module, n=integer,(default: n=0): Return (interreduced)
generators of id of degree <=1 in reduced triangular form if n=0 [non-reduced triangular form
if n!=0]
pos_def(A); A = constant, symmetric square matrix : Return int: 1 if A is positive definit , 0 if
not, -1 if unknown.
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 68
Phụ lục 4. Một số kết quả chương trình trên
cây sinh loài 4 taxa
Với mô hình Neyman thì chương trình chạy cho kết quả tốt trên những cây 4
taxa, riêng mô hình Jukes – Cantor thì chương trình chỉ chạy được với dạng cây hình
móng (hình 11).
Hình 11: Cây sinh loài 4 taxa hình móng
Với mô hình Jukes – Cantor, cây sinh loài 4 taxa hình móng (hình 11) có 5 lớp
xác suất tương đương. Lấy u = (700,7,100,100,20). Chương trình chạy cho các kết quả
29 nghiệm hợp lý, trong đó có 2 nghiệm thực dương và nghiệm làm cực đại toàn cục
hàm hợp lý là:
(0.57910415483,0.34772086852,0.036090949833,0.037084026813,0.011796766418)
và giá trị hàm hợp lý tương ứng là: -3083.02
Hình 12: Cây sinh loài 4 taxa hình cần trục
Với mô hình Neyman, cây sinh loài 4 taxa hình cần trục (hình 12) có 4 lớp xác
suất tương đương. Lấy u = (350, 7, 100, 20), chương trình chạy cho các kết quả 5
nghiệm, trong đó có 1 nghiệm thực dương:
(0.59910005385, 0.20030240782, 0.06521874649, 0.13537879184)
và giá trị hàm hợp lý tương ứng: -503.564
Sau đây một số hình dáng cây 4 taxa khác (không gốc và có gốc)
Hình 13: Một số cây sinh loài 4 taxa
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 69
Phụ lục 5. Bảng đối chiếu Thuật ngữ
Anh - Việt
Thuật ngữ Tiếng Anh Thuật ngữ Tiếng Việt
Ancestral maximum likelihood Cực đại khả năng tổ tiên
Character group Đặc trưng nhóm
Distance method Phương pháp khoảng cách
Critical point Điểm tới hạn
Evolutionary timeline Biểu đồ thời gian của tiến hóa
Fossil Hóa thạch
Homologous characters Đặc tính đồng dạng
Homogeneous polynomial Đa thức thuần nhất
Idean Iđêan
Likelihood locus Quỹ tích hợp lý
Maximum Likelihood Hợp lý cực đại
Maximum Likelihood Estimation Ước lượng hợp lý cực đại
Maximum parsimony Hà tiện tối đa
Neighbor-joining Liên kết cận kề
Paleontology Khảo cổ học
Phylogenetic invariant Bất biến trên cây sinh loài
Phenetics Ngoại hình học
Record Di chỉ
Singular point Điểm kỳ dị
Species Các loài
Tangent space Không gian tiếp xúc
Taxa Nút lá cây sinh loài
Variety Đa tạp
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 70
Danh mục các tên
Abel, 12, 13, 16, 39, 40
AML, 33
ancestral, 33
ánh xạ, 12, 37, 38, 40, 44
Bậc của đa thức, 15
bất biến, 10, 11, 36, 37, 38, 39, 42, 45,
46, 47, 48, 51, 54
bất biến của cây sinh loài, 11, 38, 39,
47
bất biến tầm thường, 46, 51
Bayes, 19
biến cố, 18, 19, 24
Biến cố chắc chắn, 18
Biến cố không thể, 18
Biến cố ngẫu nhiên, 18
biến cố sơ cấp, 18
biến ngẫu nhiên rời rạc, 20, 21
Cây sinh loài, 3, 30, 32, 36, 42
cladistics, 30
Đa tạp đại số, 18
đa tạp tuyến tính, 18
đa thức thuần nhất, 15, 46, 47, 48, 50
đặc trưng, 10, 20, 21, 22, 40, 42
đại lượng ngẫu nhiên, 10, 20
đẳng cấu nhóm, 12
dạng tuyến tính, 15
det, 15
điểm kỳ dị, 47, 50
Định thức, 15
Distance method, 31, 69
DNA, 3, 9, 30, 31, 32, 38, 39
Độc lập tuyến tính, 16
độc lập tuyến tính tối đại, 17
đồng cấu nhóm, 12, 40
evolutionary timeline, 30
fossil, 30
Fourier, 38, 40, 41, 42, 44, 46, 64
gen, 3
genomics, 30
Gradient, 47
Hà tiện tối đa, 31, 69
hàm hợp lý, 23, 26, 27, 47, 49, 51, 52,
55, 56
hàm phân phối, 10, 20
hạng của hệ vector, 17
hạng tử, 14, 15
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 71
Hạt nhân, 12
Hessian, 50, 51, 52, 66, 67
homologous characters, 9
iđêan, 13, 14, 39, 47, 48, 50
Jacobi, 47, 48
Jukes - Cantor, 34, 35
Ker, 49, 50, 62
không gian vector, 16, 17
Không gian vector, 16
Kimura, 10, 30, 35, 39, 42
Kỳ vọng, 20, 21
Lagrange, 49, 50, 51, 52, 66
likelihood, 26, 27, 28, 33, 59, 69
ln, 28
logarit, 26, 50
lớp xác suất tương đương, 43, 55
ma trận, 9, 10, 15, 17, 31, 33, 34, 35,
36, 38, 39, 43, 48, 49, 50
Ma trận chuyển vị, 15
mẫu ngẫu nhiên, 21, 22, 23, 24
Maximum parsimony, 9, 31, 69
MLE, 3, 10, 25, 26, 27, 28
neighbor-joining, 31
Neyman, 10, 30, 39
nhóm, 3, 9, 10, 11, 12, 13, 14, 16, 19,
30, 39, 40, 41, 42, 69
nhóm con, 13
nucleotide, 33, 34, 35
paleontology, 30
phân phối xác suất, 20
phenetics, 30
phép thử, 18
phép thử ngẫu nhiên, 18
phụ thuộc tuyến tính, 16, 17
phương pháp hợp lý cực đại, 9, 10, 22,
37
Phương sai, 20, 21
phương trình hợp lý, 11, 23, 24, 37, 47,
50, 51, 54, 57
phylogenetics, 30
Quỹ tích hợp lý, 47, 69
record, 30
species, 30
Tích của 2 biến cố, 18
tổ hợp tuyến tính, 17
Tổng của 2 biến cố, 18
trường, 1, 10, 14, 16, 17, 22, 26, 27, 33,
35, 36, 37, 43, 55, 56, 57
Ước lượng điểm, 22
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 72
ước lượng hợp lý cực đại, 3, 10, 12, 22,
23, 24, 25, 26, 57, 58
Ước lượng không chệch, 22
Ước lượng tham số, 22
vành con, 13
vành đa thức, 10, 14, 39, 47
vành giao hoán, 13, 14
vành thương, 14
vector, 16, 17, 47, 48, 49, 50, 54
Vector vận tốc, 47
xác định âm, 50, 51, 52
Các file đính kèm theo tài liệu này:
- 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ỏ.pdf