ĐẠI HỌC QUỐC GIA HÀ NỘI
TRƯỜNG ĐẠI HỌC CÔNG NGHỆ

Tạ Thị Thúy
Xây dựng mô hình biến đổi amino axit cho nấm
KHOÁ LUẬN TỐT NGHIỆP ĐẠI HỌC HỆ CHÍNH QUY

Ngành: Khoa học máy tính
HÀ NỘI – 2010
ĐẠI HỌC QUỐC GIA HÀ NỘI
TRƯỜNG ĐẠI HỌC CÔNG NGHỆ

Tạ Thị Thúy
Xây dựng mô hình biến đổi amino axit cho nấm
KHOÁ LUẬN TỐT NGHIỆP ĐẠI HỌC HỆ CHÍNH QUY

Ngành: Khoa học máy tính
Cán bộ hướng dẫn: TS.Lê Sỹ Vinh
Cán bộ đồng hướng dẫn: ThS.Đặng Cao Cường
HÀ NỘI - 2010
Lời cảm ơn
Đầu tiên, em xin bày tỏ lòng kính trọng và cảm ơn sâu sắc tới Tiến sỹ Lê Sĩ Vinh
và Thạc sĩ Đặng Cao Cường đã trực tiếp giao đề tài và tận tình hướng dẫn em trong
suốt quá trình thực hiện khóa luận này.
Em cũng xin bày tỏ lòng biết ơn đến các thầy cô trong trường Đại học Công nghệ
đã giảng dạy và cho em những kiến thức quý báu, làm nền tảng để hoàn thành khóa
luận cũng như thành công trong nghiên cứu, làm việc trong tương lai.
Cuối cùng, cho em gửi lời cảm ơn sâu sắc tới gia đình đã luôn động viên và tạo
điều kiện để em học tập tốt và hoàn thành được khóa luận.
Hà Nội, tháng 05, năm 2010
Sinh viên

T A C A T C G C C A T G
- Chuỗi ARN :
AAUGACUUCUAGCCGA
Thông tin di truyền được truyền từ thế hệ này sang thế hệ khác nhờ quá trình sao
chép ADN (chuỗi con sinh ra giống hệt chuỗi mẹ).
1.1.2 Amino axit và quá trình tổng hợp protein
Amino axit là đơn phân tử cấu tạo nên protein. Giống như ADN, tính đặc thù và
đa dạng của protein thể hiện bởi số lượng , thành phần và trật tự sắp xếp các amino
axit trong chuỗi. Có tất cả 20 loại amino axit cấu tạo nên các protein.
Chuỗi amino axit được tổng hợp từ một đoạn của chuỗi ADN (gen) thông qua
quá trình sau:
ADN ARN protein
Trong quá trình sao mã (transcription), trình tự các nucleotid trong một chuỗi
đơn của chuỗi ADN, gọi là mạch gốc, sẽ quy định trình tự các ribonucleotid trong
6
Dịch mã
Sao mã
ARN theo nguyên tắc bổ sung (A = U, G = C). Trong quá trình dịch mã (translation),
cứ bộ ba ribonucleotid trong ARN sẽ quy định 1 amino axit trong protein. Hay nói
cách khác, mỗi bộ ba nucleotid trong chuỗi ADN (được gọi là một codon hay mã di
truyền) sẽ mã hóa cho 1 amino axit trong chuỗi protein được tổng hợp.
Ví dụ về quá trình sinh tổng hợp protein:
- Mạch gốc ADN: T A C G C C A A G A T T
Sao mã
- Mạch ARN A U G C G G U U C U A A
Dịch mã
- Chuỗi amino axit: aa
mđ
aa
1

Glycine Gly G GGU, GGC, GGA, GGG
Histidine His H CAU, CAC
Isoleucine Ile I AUU, AUC, AUA
Lysine Lys K AAA, AAG
Leucine Leu L UUA, UUG, CUU, CUC, CUA, CUG
Methionine Met M AUG (START)
Asparagine Asn N AAU, AAC
Proline Pro P CCU, CCC, CCA, CCG
Glutamine Gln Q CAA, CAG
Arginine Arg R CGU, CGC, CGA, CGG, AGA, AGG
Serine Ser S UCU, UCC, UCA, UCG, AGU, AGC
Threonine Thr T ACU, ACC, ACA, ACG
Valine Val V GUU, GUC, GUA, GUG
Trytophan Trp W UGG
Tyrosine Tyr Y UAU, UAC
STOP UAA, UGA, UAG
Bảng 1. Mã univesal của 20 loại amino axit
1.1.3 Đột biến và tương đồng
Trong quá trình phát triển và tiến hóa, dưới tác động của môi trường và các tác
nhân hóa học, các quá trình đột biến xảy ra, trải qua chọn lọc tự nhiên và được tích lũy
dần dần dẫn đến sự biến đổi thông tin di truyền giữa các thế hệ, hình thành nên sự đa
dạng sinh học. Từ chuỗi ADN ban đầu sau các quá trình đột biến sẽ tạo ra các chuỗi
ADN con không còn giống hệt chuỗi ADN mẹ nữa, nếu được chọn lọc tự nhiên chấp
nhận nó sẽ tồn tại và truyền cho đời sau như một chuỗi mới độc lập. Sự kiện đó gọi là
sự phân kỳ chuỗi. Các đột biến được xét ở đây là các đột biến điểm bao gồm:
- mất 1 nucleotid.
- thêm 1 nucleotid.
- thay thế nucleotid này thành nucleotid khác.
Việc biến đổi chuỗi ADN (tức là biến đổi bộ mã hóa cho chuỗi protein) có thể
dẫn đến làm biến đổi chuỗi amino axit tương ứng:

1
0
1
1
1
2
D
1
E I H K I R M T L T S T
D
2
E I L K I R I T L - S T
D
3
- I L K I R I T L T S T
D
4
V I H K I R I N L T S T
D
5
- I L K I R M T L - S T
Trong ví dụ trên, alignment chỉ ra đột biến thay thế amino axit xảy ra tại vị trí
số 1,3,7 và 8; và đột biến thêm hoặc mất amino axit xảy ra ở vị trí số 1 và 10.
1.1.5 Khoảng cách tiến hóa
Khi có hai chuỗi amino axit (nucleotid) tương đồng, vấn đề đặt ra làm thế nào để
biết được mối quan hệ “họ hàng” giữa chúng xa hay gần như thế nào. Độ đo mối quan
hệ này gọi là khoảng cách tiến hóa (evolutionary distance hay genetic distance), nó
biểu thị số biến đổi xảy ra giữa hai chuỗi.
Cách đơn giản nhất để đo khoảng cách tiến hóa giữa hai chuỗi tương đồng là sắp
hàng cho chúng và đếm số vị trí tương đồng khác biệt. Tỉ lệ khác biệt này còn được

dùng khoảng cách quan sát để thể hiện khoảng cách tiến hóa không còn chính xác nữa.
Nhìn chung, khoảng cách quan sát sẽ nhỏ hơn khoảng cách tiến hóa. Có thể suy ra
khoảng cách tiến hóa từ khoảng cách quan sát nhưng mối liên quan này thường không
tuyến tính.
Khoảng cách tiến hóa cũng có thể suy ra trực tiếp từ dữ liệu bằng phương pháp
maximum-likelihood (sẽ được trình bày ở phần sau).
1.1.6 Cây phát sinh loài
Nếu khoảng cách tiến hóa chỉ thể hiện mối quan hệ giữa hai chuỗi tương đồng,
thì cây phát sinh loài có thể minh họa mối quan hệ tiến hóa đồng thời giữa tất cả các
chuỗi (hay giữa các loài) trong một tập dữ liệu. Trong đó, lá của cây là các loài trong
thời kì hiện tại (hay là các chuỗi DNA hoặc protein tương ứng đại diện cho mỗi loài).
Các nút trong được coi như là các tổ tiên giả định. Các nhánh chính là khoảng cách
tiến hóa giữa hai loài (hay nói cách khác là giữa hai chuỗi). Cây phát sinh loài thường
là cây có các nhánh rẽ đôi với giả thuyết là từ một chuỗi tổ tiên ban đầu chỉ có thể làm
xuất hiện 2 chuỗi con mới phát triển độc lập. Cây tiến hóa có thể có gốc hoặc không có
gốc. Nếu là cây có gốc, khi đó gốc cây sẽ được coi là tổ tiên chung của tất cả các loài
trong cây. Nếu cây không có gốc tức là cây chỉ thể hiện mối quan hệ gần gũi giữa các
loài trong cây.
Hình 1 là một ví dụ đơn giản về một cây phát sinh loài trong đó các lá A, B, C,
D là các loài hiện tại, các nút trong là các tổ tiên giả định, các số là độ dài nhánh cây.

1
Hình 1. Ví dụ về cây phát sinh loài
Việc xây dựng cây phát sinh loài đều cần đến việc sắp hàng trình tự các chuỗi
và một mô hình tiến hóa.
12
2
B
A
1

ii
được đặt giá trị sao cho tổng
mỗi hàng bằng 0.
Mô hình thừa nhận 4 giả định:
- Tốc độ thay thế từ nucleotid i thành nucleotid j là độc lập với các quá trình
trước đó sinh ra nucleotid i (tính chất Markov)
- Tốc độ biến đổi là không đổi theo thời gian (time-homogeneous).
- Đột biến thay thế giữa các nucleotid có thể xảy ra ở bất kỳ thời gian nào của
quá trình (time-continuous).
- Tần suất π = (π
A
, π
C
, π
G
, π
T
) của các nucleotid là cân bằng (stationary).
Một số mô hình còn có thêm giả định time-reversibility, tức là tốc độ biến đổi từ
nucleotid i thành nucleotid j bằng tốc độ ngược lại biến đổi từ nucleotid j về nucleotid
i. Hay tốc độ tương đối sẽ là a'=a, b'=b, c'=c, d'=d, e'=e và f'=f. Mô hình như vậy gọi là
GTR (general time-reversible model). Khi đó mô hình Q trở thành:
A C G T

13
Ma trận Q có thể phân tích thành ma trận tốc độ biến đổi tương đối hay ma trận
tốc độ trao đổi giữa các amino axit R = {R
ij
} và tần suất nucleotid π như sau:
Trong đó ma trận tốc độ trao đổi có dạng:

axit chứa nhiều tham số và để ước lượng được mô hình thì khối lượng tính toán rất
lớn. Do đó các mô hình biến đổi amino axit thường là các mô hình kinh nghiệm tức là
các mô hình được ước lượng một lần trên tập dữ liệu lớn và được sử dụng lại vào các
bài toán cụ thể.
1.2.3 Ước lượng khoảng cách tiến hóa
Khi đã có mô hình biến đổi Q, tức là các tham số về tần suất và tốc độ biến đổi
của từng loại amino axit đã biết, thì có thể tính được khoảng cách tiến hóa giữa 2
chuỗi amino axit bằng phương pháp maximum likelihood. Hàm likelihood L(d) là hàm
tính xác suất biến đổi từ 1 chuỗi thành chuỗi còn lại nếu trung bình có d biến đổi xảy
ra trên mỗi vị trí (hay khoảng cách tiến hóa là d), cụ thể:
Trong đó, x
i
là amino axit tại vị trí thứ i trong chuỗi x = {x
i
}. P
xiyi
(d) là xác suất
chuyển từ amino axit x
i
trong chuỗi x sang amino axit y
i
trong chuỗi y sau thời gian
tiến hóa d. Khi đó d được ước lượng bằng cách tìm giá trị d* làm cực đại hóa hàm
likelihood L(d):
d* = argmax {L(d)}
1.2.4 Mô hình với tốc độ biến đổi theo vị trí
Thực nghiệm cho thấy tốc độ xảy ra đột biến thay thế là không giống nhau giữa
các vị trí. Ví dụ, tốc độ biến đổi ở vị trí thứ ba của 1 codon trong chuỗi nucleotid mã
hóa cho protein thường nhanh hơn ở vị trí thứ nhất và thứ hai. Đối với chuỗi protein,
các vị trí mà nếu biến đổi nó sẽ ít tác động lên chức năng hoặc cấu trúc protein thì tốc

số α nhỏ hơn mô tả sự biến đổi tốc độ theo vị trí lớn hơn.
16
Hình 2. Phân phối Gamma
Hàm Γ(α) liên tục có thể được thay xấp xỉ bằng hàm Γ rời rạc với c lớp tương
ứng với c hệ số tỉ lệ tốc độ biến đổi r
1
, r
2
, , r
c
(Yang, 1994). Tham số α thường được
ước lượng từ dữ liệu.
Có thể kết hợp mô hình 2 trạng thái với mô hình sử dụng phân phối Γ. Khi đó,
mô hình lai sẽ giả định tỉ lệ θ các vị trí là không đổi, các vị trí khác là biến đổi với tốc
độ biến đổi tuân theo phân phối Γ.
1.2.5 Tại sao sử dụng mô hình biến đổi amino axit
Mô hình biến đổi amino axit là một trong những thành phần quan trọng trong các
bài toán phân tích sinh học liên quan đến protein như:
- Sắp hàng trình tự chuỗi protein
- Tính khoảng cách tiến hóa giữa các chuỗi
- Xây dựng cây phát sinh loài (cây tiến hóa).
- Các ứng dụng khác như dự đoán chức năng của protein mới,
Đối với bài toán xây dựng cây phát sinh loài, tùy dữ liệu cần phân tích là các
chuỗi nucleotid hay amino axit, mô hình tiến hóa được sử dụng là mô hình biến đổi
nucleotid hay mô hình biến đổi amino axit. Việc suy ra cây tiến hóa từ chuỗi amino
axit có một vài ưu điểm so với chuỗi nucleotid [10] như:
Thứ nhất, gen được thể hiện chỉ khi nó là bộ mã hóa cho chuỗi protein và ở dạng
tập các bộ ba mã hóa (gọi là gen cấu trúc). Ở các động vật bậc cao, các gen có thể
chứa các đoạn intron, tức các đoạn không mã hóa protein nằm xen kẽ giữa các exon –
17

Maximum likelihood
MtMam [4] 1998 Protein mã hóa bởi chuỗi DNA
nằm trong ti thể
Maximum likelihood
CpREV [3] 2000 Protein mã hóa bởi chuỗi DNA
nằm trong lục lạp *** (thực
vật)
Maximum likelihood
VT [12] 2000 Mô hình chung resolvent
WAG [15] 2001 Mô hình chung Maximum likelihood
rtREV [6] 2002 Retroviral Pol protein Maximum likelihood
DcMut [9] 2004 Mô hình chung Phương pháp đếm
MtArt [1] 2006 Protein mã hóa bởi chuỗi DNA
nằm trong ti thể (động vật thân
đốt)
Maximum likelihood
HIVb [13] 2007 Vi rút HIV Maximum likelihood
HIVw [13] 2007 Vi rút HIV Maximum likelihood
19
LG [14] 2008 Pfam – mô hình chung Maximum likelihood
FLU [5] 2009 Vi rút cúm Maximum likelihood
Bảng 2. 15 mô hình biến đổi amino axit hiện tại
(*) Mô hình chung có nghĩa là tập dữ liệu lớn và phân kỳ gồm các chuỗi protein
thuộc nhiều loài, nhiều họ protein khác nhau.
(**) ti thể là một cơ quan trong tế bào.
(***) lục lạp là một cơ quan trong tế bào thực vật.
2.2. Mô hình Dayhoff, JTT
Dayhoff và các cộng sự [11] là những người đầu tiên mô hình hóa quá trình biến
đổi amino axit (1978). Họ sử dụng tập dữ liệu gồm 71 tập các protein có họ hàng gần
nhau và quan sát được 1572 biến đổi giữa các amino axit. Bằng phương pháp đếm đơn

Năm 2001, Whelan và Goldman [15] cải tiến phương pháp maximum likelihood
để ước lượng ra mô hình biến đổi amino axit mới, gọi là mô hình WAG. Hai ông chỉ ra
rằng có thể sử dụng các cây tiến hóa gần đúng để ước lượng ra mô hình biến đổi với
kết quả không ảnh hưởng nhiều, do đó có thể đơn giản hóa việc tính toán bằng cách
tránh việc ước lượng đồng thời tham số mô hình, cây tiến hóa và độ dài nhánh cây. Họ
sử dụng tập dữ liệu lớn hơn, BRKALN, gồm 182 họ protein và ~900000 amino axit.
Các bước ước lượng cây bao gồm: trước hết, các cây phát sinh loài được khởi tạo bằng
phương pháp NJ (Neighbor-Joining) – một phương pháp xây dựng cây nhanh, đơn
giản; sau đó, ước lượng lại các độ dài nhánh cây sử dụng phương pháp maximum
likelihood với mô hình khởi tạo là JTT; tiếp tục, sử dụng cây này để ước lượng tham
số mô hình mới bằng phương pháp maximum likelihood. Mô hình WAG cho thấy tốt
hơn mô hình JTT và mô hình Dayhoff dựa trên so sánh giá trị likelihood của các cây
phát sinh loài, tức là các cây phát sinh loài được khởi tạo sử dụng mô hình WAG mô
phỏng mối quan hệ tiến hóa của dữ liệu tốt hơn.
2.6. Mô hình LG
Mô hình LG (Le & Gascuel, 2008) hiện nay được coi là mô hình chung tốt nhất.
Cũng cách tiếp cận maximum likelihood như mô hình WAG nhưng Le và Gascuel đã
cải tiến bằng cách xét thêm yếu tố biến đổi tốc độ theo vị trí trong quá trình ước lượng
mô hình và sử dụng tập dữ liệu lớn hơn, Pfam (gồm ~50000 chuỗi, ~ 6,5 triệu amino
axit).
21
2.7. Mô hình FLU
Các năm gần đây chứng kiến sự xuất hiện và phát triển nhiều dịch bệnh mà
nguyên nhân là virut, do đó các bài toán phân tích quá trình phát triển và tiến hóa của
virut yêu cầu độ chính xác cao hơn. Cũng áp dụng cách tiếp cận maximum likelihood,
C.Dang, Q.Le, O.Gascuel và V.Le đã ước lượng ra mô hình mới (FLU) cho loài virut
cúm [5]. Dữ liệu sử dụng gồm ~113,000 chuỗi protein thuộc các chủng virut cúm với
tổng số ~20 triệu amino axit. Mô hình FLU chứng tỏ nó tốt hơn hẳn 14 mô hình trước
đó khi mô phỏng quá trình tiến hóa của virut cúm.
22

=
{q
xy
} được giữ không đổi trong suốt quá trình tiến hóa. Trong đó q
xy
(x#y) là số biến
đổi từ amino axit x thành amino axit y xảy ra trong một đơn vị thời gian, phần tử
đường chéo q
xx
được gán giá trị sao cho tổng mỗi hàng bằng 0. Mô hình thừa nhận các
23
giả định: time-reversible, time-continuous, time-homogeneous và stationary (xem
chương 1). Khi đó, ma trận Q có thể phân tích thành ma trận tốc độ trao đổi R = {r
xy
}
và tần suất amino axit π = {π
xy
}, trong đó:
q
xy
= π
y
r
x↔y
với x # y
q
xx
= -∑ q
xy
với x = y (1)

i
Trong đó L(T,Q|D
i
) là xác suất của dữ liệu tại ví trí i trong alignment D được cho
bởi cây T và mô hình Q. L(T,Q|D
i
) được tính bằng cách áp dụng phương trình (2) cho
mỗi nhánh cây (t là độ dài nhánh) và sử dụng giải thuật prunning [10].
Phương trình (3) tính xác suất của D với giả định rằng tốc độ biến đổi là như
nhau giữa các vị trí. Đây là giả định không thực tế và để khắc phục nó, các mô hình
khác nhau mô hình hóa sự biến đổi tốc độ theo vị trí đã được đề xuất, trong đó hay
được sử dụng là mô hình hai trạng thái và mô hình phân phối Γ (xem chương 1). Mô
hình kết hợp của hai mô hình trên (kí hiệu mô hình r) giả định một giá trị θ
inv
là tỉ lệ số
24
vị trí trong chuỗi không biến đổi, các vị trí khác biến đổi với tốc độ biến đổi phân bố
theo phân phối Γ. Khi đó xác suất hay hàm likelihood của dữ liệu khi tuân theo cây
tiến hóa T, mô hình tiến hóa Q và mô hình biến đổi tốc độ theo vị trí r là:
(4)
Trong đó L(invariant|d
i
) là xác suất của vị trí d
i
theo mô hình hai trạng thái. Theo
đó, L(invariant|d
i
) bằng π
x
nếu vị trí d

gần với giá trị tối ưu. Khi đó việc tính toán được đơn giản hóa và Q được ước lượng
xấp xỉ bằng cách cực đại hóa hàm likelihood:
(6)
Quá trình ước lượng mô hình được thể hiện dưới sơ đồ sau:
25