EM 알고리즘을 통한 haplotype frequency 추정


우선 EM 알고리즘의 정의부터 알아보자.


EM 알고리즘은 관측되지 않는 잠재변수(unobserved latent variables)에 의존하는 확률 모델에서 최대우도(maximum likelihood)최대사후확률(Maximum A Posteriori)을 갖는 매개변수를 찾는 반복적인 알고리즘이다.


[출처] EM 알고리즘(1/2)|작성자 lhm0812


그렇다면 왜 haplotype frequency 추정할 때 EM 알고리즘 쓸 수 있는가? 관측되지 않은 변수가 있기 때문이다. 모수를 추정하고자 하는데 모수 추정에 필요한 값이 관측할 수 없는 잠재변수를 포함하고 있을 때 쓸 수 있다.


아래와 같은 세팅을 생각해보자.


Locus 1 has alleles A; a
Locus 2 has alleles B; b


Genotype counts are observed at the two loci for N individuals in the population.


homozygous for both loci: nAABB; naaBB ; nAAbb; naabb
homozygous for one loci: nAABb; naaBb; nAaBB; nAabb
heterozygous for both loci: nAaBb


근데 문제는 여기서 nAaBb는 2개의 가능한 하플로타입 nAB/ab or nAb/aB 이 있다. (즉, 두 loci 모두 genotype이 heterozygous인 경우 haplotype을 결정할 수 없다는 사실을 알 수 있다.)


우리가 알고 싶은건 pAB; pAb; paB; pab. 바로 이 haplotype frequency들이다.

여기서 pAB를 예로 들어 EM 알고리즘을 적용하는 과정을 살펴보자. 우선 pAB의 MLE는 아래와 같다.



p^AB = nAB/2N (N=관측된 사람 수. 하플로타입은 개인당 2개이므로 2N)



또한 nAB = 2nAB/AB + nAB/Ab + naB/AB + nAB/ab


이 때 문제는 nAB/ab를 모른다. 이를 잠재변수라고 하고, 이는 추정하려는 모수들 (pAB; pAb; paB; pab)과 관련되어 있다. 그래서 모수가 있다면 구할 수 있는 값이다. EM 알고리즘에서는 이를 모수를 적절한 값으로 초기화한 후, 반복적으로 nAB/ab를 구한다.



EM 알고리즘에서는 추정하려는 모수 pAB; pAb; paB; pab 를 적절한 값(ex. 0.25)으로 초기화한다. 그리고 이를 p0AB 처럼 표기한다.


관측한 값

Y = (nAABB; naaBB; nAAbb; naabb; nAABb; naaBb; nAaBB; nAabb; nAaBb)


와 초기화된 모수 p가 주어졌을 때, nAB의 기댓값은 아래와 같이 구할 수 있다.




여기서 E[nAB/ab]는 무엇인가?



왜 이렇게 되냐하면, nAB/ab는 nAaBb는 nAB/ab 또는 nAb/aB 이므로, 두 하플로타입(AB/ab, Ab/aB)이 등장하는 확률 값의 비율을 통해 nAB/ab의 기댓값을 구할 수 있다.


이를 통해 n0AB를 구했으면 pAB의 MLE를 통해 이를 업데이트한다.


p^AB = nAB/2N


이 과정을 반복하면 p^AB가 수렴해가고 이를 통해 추정값을 얻을 수 있다. 다른 모수들에 대해서도 이와 같이 한다.



출처 - http://courses.washington.edu/b516/lectures_2010/EM_Algorithm_Haplotype_Frequency_2010.pdf


'Domains > Genetics' 카테고리의 다른 글

Genetics - Genetic Mapping1  (0) 2017.11.11
Genetics - Recombination  (0) 2017.11.11
EM 알고리즘 haplotype frequency 추정  (0) 2017.10.08
DNA 결합 단백질  (2) 2017.08.12
유전체학을 위한 신경망 모델  (0) 2017.08.12
Regulatory Genomics 개념과 ENCODE 프로젝트  (0) 2017.08.08