BioPlayground

🧬
목록으로

계통수 최대우도법: IQ-TREE와 RAxML의 통계 모델·수치 최적화

파시모니의 함정을 넘어선 통계적 계통수. JC69부터 GTR+Γ+I까지 진화 모델, Felsenstein의 pruning 알고리즘, IQ-TREE의 실전 사용까지.

심화
|
20분
|
검증 완료 (2026-07-20)
phylogenymaximum likelihoodsubstitution model
진행률0/34 (0%)

파시모니의 한계에서 시작하기

M27의 마지막에 파시모니의 함정 "long branch attraction"을 언급했습니다. 최대우도(Maximum Likelihood, ML)는 이 함정을 통계적으로 넘어섭니다. 진화 모델을 명시적으로 정의하고, 그 모델 하에서 관측 데이터의 확률이 가장 큰 트리를 찾습니다.

이 편에서는 진화 치환 모델의 계층, Felsenstein의 pruning 알고리즘, IQ-TREE 실전 사용, 그리고 SH-aLRT·UFBoot 같은 현대 트리 신뢰도 지표까지 훑습니다. Micro Tier 계통수 3편의 왕관이자 실전 표준입니다.

진화 치환 모델 계층

DNA 서열의 치환은 확률 과정으로 표현됩니다. Rate matrix Q로 정의.

JC69 (Jukes-Cantor): 모든 치환률 동일. 파라미터 1개.

Q=(−αααα−αααα−αααα−)Q = \begin{pmatrix} - & \alpha & \alpha & \alpha \\ \alpha & - & \alpha & \alpha \\ \alpha & \alpha & - & \alpha \\ \alpha & \alpha & \alpha & - \end{pmatrix}

K80 (Kimura): transition과 transversion 구분. 파라미터 2개.

HKY85: transition/transversion + 4자 빈도 불균등. 파라미터 5개.

GTR (General Time Reversible): 6가지 치환률 + 4자 빈도. 파라미터 9개. 가장 유연.

GTR+Γ: GTR + 자리별 rate 이질성(감마 분포). 실무 default.

GTR+Γ+I: + invariant sites 비율. 극도로 이질적 데이터에.

IQ-TREE의 -m TEST 옵션이 자동으로 최적 모델을 선택합니다. 실무에서 모델 선택을 손으로 하는 시대는 지났습니다.

트리 우도 = 관측 MSA의 확률

트리 T, 브랜치 길이 b, 진화 모델 M이 주어졌을 때, 관측 MSA X의 우도는

L(T,b,M)=P(X∣T,b,M)=∏site iP(columni∣T,b,M)L(T, b, M) = P(X | T, b, M) = \prod_{\text{site } i} P(\text{column}_i | T, b, M)

각 컬럼은 독립이라 가정. 로그를 취하면

log⁡L=∑ilog⁡P(columni∣T,b,M)\log L = \sum_i \log P(\text{column}_i | T, b, M)

한 컬럼의 확률 P(column | T, b, M) 계산이 문제입니다. 각 내부 노드에 상태가 미상이므로, 모든 가능한 조상 상태 조합을 합해야 합니다. 정공법은 상태 수의 (내부 노드 수)제곱 = 지수 시간. 감당 불가.

Felsenstein의 pruning 알고리즘이 이 문제를 다항 시간으로 낮췄습니다.

Felsenstein Pruning — HMM의 Forward-Backward와 같은 뼈대

M19에서 HMM Forward를 배웠습니다. 트리 위의 pruning은 그 확장입니다.

각 노드 v와 각 상태 x에 대해 조건부 우도 L_v(x) 정의:

Lv(x)=P(v 아래 서브트리의 관측∣v의 상태=x)L_v(x) = P(v \text{ 아래 서브트리의 관측} | v \text{의 상태} = x)

Leaf: 관측 문자가 x면 L_leaf(x) = 1, 아니면 0.

Internal node (자식 c1, c2, 브랜치 길이 b1, b2):

Lv(x)=[∑yP(x→y∣b1)Lc1(y)]×[∑yP(x→y∣b2)Lc2(y)]L_v(x) = \left[ \sum_{y} P(x \to y | b_1) L_{c_1}(y) \right] \times \left[ \sum_{y} P(x \to y | b_2) L_{c_2}(y) \right]

P(x → y | b)는 브랜치 길이 b 동안 상태 x가 y로 변할 확률. 진화 모델 M에서 파생 (Q의 지수).

Root에서:

P(column)=∑xπ(x)Lroot(x)P(\text{column}) = \sum_x \pi(x) L_{\text{root}}(x)

π(x)는 초기 상태 확률(equilibrium frequency).

이 계산은 각 컬럼에 대해 O(N × K²). N은 트리 크기, K는 상태 수(DNA=4, 단백질=20). 각 컬럼을 독립 병렬 처리하면 대규모 MSA도 실전 가능. HMM Forward가 시간 축을 따라 확률을 축적한다면, Felsenstein pruning은 트리 축을 따라 확률을 축적합니다. 같은 DP 뼈대.

브랜치 길이 최적화

트리 위상 고정 상태에서 브랜치 길이 b를 최적화합니다. 로그 우도를 각 b로 미분하고 뉴턴법이나 golden section search로 수렴. 각 브랜치 최적화 후 다시 다음 브랜치 최적화. 여러 iteration 반복. 이 부분이 실전 도구의 수치 최적화 핵심입니다.

트리 위상 탐색

브랜치 길이는 다항 시간이지만, 트리 위상 자체는 large ML도 파시모니처럼 NP-완전. 실전 도구는 heuristic search를 씁니다.

  • NNI (Nearest Neighbor Interchange): 이웃 브랜치 스왑
  • SPR (Subtree Pruning and Regrafting): 서브트리 잘라내서 다른 위치에 재부착
  • TBR (Tree Bisection and Reconnection): 트리 두 조각으로 자르고 재연결

IQ-TREE의 특징은 여러 초기 트리 병렬 탐색 + tree perturbation. Local optimum 함정을 회피합니다.

IQ-TREE 실행 — 실전 표준

bash
# 설치
wget https://github.com/iqtree/iqtree3/releases/download/v3.0.0/iqtree-3.0.0-Linux.tar.gz
tar xzf iqtree-3.0.0-Linux.tar.gz
./iqtree3-3.0.0-Linux/bin/iqtree3 --version


# 기본 실행 (자동 모델 선택 + ML 트리 + 부트스트랩)
iqtree3 -s my_msa.fasta -m TEST -B 1000 -T AUTO


# 옵션 해설
# -s: 입력 MSA
# -m TEST: 자동 진화 모델 선택 (또는 -m GTR+G 직접 지정)
# -B 1000: UFBoot 부트스트랩 1000회
# -T AUTO: CPU 스레드 자동

주요 출력 파일:

  • *.treefile: Newick 포맷 최대우도 트리
  • *.iqtree: 실행 리포트 (선택된 모델 · 우도 · 부트스트랩 요약)
  • *.log: 상세 로그

현대 신뢰도 지표 — UFBoot과 SH-aLRT

M26에서 표준 부트스트랩을 소개했지만, IQ-TREE는 두 가지 최신 지표를 씁니다.

UFBoot (Ultrafast Bootstrap): Minh et al. 2013. RELL(RandomELL) 재추정으로 표준 부트스트랩 대비 10～100배 빠름. 값 해석은 표준 부트스트랩과 살짝 다름: UFBoot ≥ 95이 신뢰 임계 (표준 부트스트랩의 70에 해당).

SH-aLRT (Shimodaira-Hasegawa approximate Likelihood Ratio Test): 각 브랜치가 통계적으로 유의한지 우도비 검정. SH-aLRT ≥ 80이 유의 임계.

두 지표 모두 통과한 브랜치가 "잘 지지된" 브랜치. 실무에서 트리 그림을 그릴 때 두 값을 함께 표시합니다.

부분 시각화 예제

트리 파일이 있으면 iTOL(https://itol.embl.de/) 이나 FigTree로 시각화. 또는 파이썬 ete3 라이브러리.

python
# ete3로 트리 로딩 + 부트스트랩 표시
from ete3 import Tree


t = Tree("my_tree.treefile", format=1)
for node in t.traverse():
    if not node.is_leaf():
        # UFBoot 값이 support attribute에
        node.name = f"{node.support:.0f}"
print(t.get_ascii(show_internal=True))

왜 파시모니 대신 최대우도인가

세 가지 결정적 이점.

  1. Long branch attraction 회피: 진화 모델이 우연 매칭을 확률적으로 감가
  2. 통계적 검정 가능: 두 트리의 우도비를 검정 (SH test, AU test)
  3. 분자 시계 검정: 시계 가정을 명시적으로 넣거나 뺄 수 있음

단점은 계산 비용. 파시모니가 O(N × K)라면 최대우도는 O(N × K² × 브랜치 길이 iteration). 실무 서열 수백 개는 노트북 몇 분. 수천～수만은 서버.

복잡도

  • 한 컬럼 우도: O(N × K²)
  • 전체 MSA 우도: O(N × K² × L). L은 컬럼 수
  • 브랜치 길이 최적화: 각 브랜치 몇 iteration
  • 트리 위상 탐색: NNI/SPR 각 이동마다 우도 재계산

IQ-TREE의 병렬화가 뛰어나서 32-core로 서열 1000개 × 컬럼 500이 시간 단위.

CS 매핑 — 트리 DP의 확률판

M19의 HMM Forward, M27의 파시모니, 그리고 이 편의 Felsenstein pruning. 세 개가 정확히 같은 뼈대입니다.

  • HMM Forward: 시간 축 위의 확률 축적
  • 파시모니: 트리 위의 최소 비용 축적 (min)
  • Felsenstein pruning: 트리 위의 확률 축적 (sum)

뼈대 하나가 여러 옷을 갈아입습니다. Min을 Sum으로, 시간을 트리로 바꾸는 것만으로 알고리즘의 도메인이 완전히 달라집니다. 이 통찰이 몸에 붙으면 새 도메인의 알고리즘도 익숙한 뼈대로 즉시 해체됩니다.

Rosalind에서 채점받기

Rosalind는 최대우도 트리 문제가 직접 없지만, BA10G (Baum-Welch)를 다시 보면 Felsenstein pruning과 같은 뼈대의 다른 응용임을 확인할 수 있습니다.

다음 편으로 이어지는 갈래

  • 다음 편 (M29): Coalescent theory — 계통수를 시간 역산으로 유도. 집단 유전학 도입
  • 마지막 편 (M30): Wright-Fisher model — 세대별 확률 과정. 중성 진화의 통계 모델
  • 확장: BEAST/BEAST2 — 베이지안 계통수. 사전 확률 · MCMC 샘플링. IQ-TREE와 상보적

더 깊게 파고 싶다면

  • Felsenstein (1981), Evolutionary trees from DNA sequences: A maximum likelihood approach, J Mol Evol — Pruning 알고리즘 원 논문.
  • Minh, Nguyen, von Haeseler (2013), Ultrafast approximation for phylogenetic bootstrap, MBE — UFBoot 원 논문.
  • Nguyen, Schmidt, von Haeseler, Minh (2015), IQ-TREE: A fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies, MBE — IQ-TREE 원 논문.
  • Felsenstein Inferring Phylogenies — 표준 교재. Chapter 16～18이 ML 핵심.
  • IQ-TREE 튜토리얼: http://www.iqtree.org/doc/ — 실전 옵션 전량.
  • UC Berkeley BIDS — Yun S. Song 교수의 최대우도 강의.

IQ-TREE로 자기 관심 유전자의 상동체 30～100개를 뽑아서 트리를 그려봅시다. 파시모니 트리(MEGA로)와 위상이 어떻게 다른지 비교하면 이 편의 감이 완성됩니다.