F27～F31을 마치고도 남는 질문 — 왜 어떤 변이는 계속 놓칠까

F27～F31에서 파이프라인 순서(Snakemake·Nextflow)와 실행 환경(Docker·Slurm)을 봉인했습니다. 그런데 이 모든 인프라가 완벽해도, 정렬 단계에서 쓰는 참조 유전체 자체가 한 사람의 서열이라면 이야기가 달라집니다. GRCh38은 사실상 소수 기증자의 서열을 짜깁기한 선형(linear) 서열입니다. 이 참조와 크게 다른 대립유전자를 가진 사람의 리드는 정렬 단계에서부터 밀려나기 쉽고, 그 결과 변이 자체가 콜링되지 않고 사라집니다. F32에서는 이 문제를 그래프 구조로 풀어내는 Pangenome(범유전체) 접근과, 이를 변이 검출에 적용한 PanVariants를 다룹니다.

원리 — 참조 편향과 그래프 정렬

참조 편향이 생기는 자리

BWA-MEM2 같은 선형 정렬기는 리드를 참조 서열 하나에 맞춰 seed-and-extend 방식으로 정렬합니다. 리드와 참조 사이 불일치(mismatch)가 많아질수록 정렬 점수가 낮아지고, 어느 지점을 넘으면 정렬기는 그 위치를 아예 후보에서 제외하거나 매핑 품질(MAPQ)을 낮춰버립니다. 문제는 이 불일치가 오류가 아니라 실제 존재하는 대립유전자일 때입니다. HLA 영역처럼 인구집단 간 다형성이 큰 구간, 또는 특정 집단에서만 흔한 삽입·결실(indel)을 가진 구간에서는 참조와 다른 서열을 가진 사람일수록 정렬 실패율이 체계적으로 높아집니다. 즉 참조 편향은 무작위 오류가 아니라, 참조 서열을 만든 기증자 집단과 얼마나 가까운가에 따라 방향이 정해지는 편향입니다.

손 계산 예제: 대립유전자 빈도가 정렬 성공률에 미치는 영향

리드 길이 150bp에서 seed-and-extend 정렬기가 허용하는 최대 불일치를 6개라고 가정합시다. 참조에 없는 어떤 삽입 변이가 8bp라면, 이 삽입을 포함한 리드는 참조와의 갭 정렬 점수가 임계값 아래로 떨어져 정렬 자체가 실패하거나 다른 위치로 잘못 정렬(mismap)됩니다. 이 변이의 인구집단 내 대립유전자 빈도가 40%라면, 이형접합자 기준으로 그 위치를 커버하는 리드의 절반가량이 이 실패 경로를 탑니다.

$P(\text{정렬 실패} \mid \text{변이 보유}) \approx 1,\quad P(\text{변이 보유}) = 2pq \;(p=0.4,\ q=0.6) \approx 0.48$

이형접합자 절반에서 이 지점의 커버리지가 사실상 반토막 나고, GATK4의 변이 콜러가 최소 커버리지 임계값을 통과하지 못하면 변이 자체가 VCF에서 사라집니다. 참조 서열이 그 대립유전자를 아예 모른다는 사실 하나가, 통계적 검정력이 아니라 정렬 단계에서 변이를 지워버리는 셈입니다.

Pangenome 그래프 — 대안 서열을 그래프의 경로로 흡수하기

Pangenome(범유전체) 그래프는 참조를 서열 하나가 아니라, 알려진 변이를 반영한 그래프로 표현합니다. 노드는 서열 조각이고, 간선은 그 조각들이 실제로 이어질 수 있는 경로를 뜻합니다. SNP 하나는 노드가 갈라졌다 다시 합쳐지는 작은 버블(bubble)로, 삽입·결실은 우회 경로로, 구조 변이는 더 복잡한 서브그래프로 표현됩니다. Human Pangenome Reference Consortium이 2023년 Nature에 발표한 초안 범유전체는 47명(94 하플로타입)의 완전한 조립 서열을 하나의 그래프로 통합해, 기존 GRCh38이 대표하지 못하던 구조 변이 상당수를 그래프 안에 명시적으로 담았습니다.

vg giraffe와 GBWT — 그래프 위에서 빠르게 정렬하기

그래프가 커지면 그래프 위에서 정렬하는 비용도 커집니다. vg giraffe는 이 문제를 두 가지로 완화합니다. 첫째, 미니마이저(minimizer) 기반 시딩으로 후보 위치를 그래프 전체가 아니라 좁은 영역으로 좁힙니다. 둘째, GBWT(Graph Burrows-Wheeler Transform)로 그래프 안에 포함된 개별 하플로타입 경로를 압축 인덱싱해, 실제로 관찰된 적 없는 재조합 경로(예: 서로 다른 하플로타입의 조각이 뒤섞인 경로)로는 정렬이 새지 않도록 유도합니다. GBWT는 개념적으로 M14～M15에서 다룬 BWT·FM-index를 그래프·다중 서열로 확장한 자료구조입니다.

PanVariants — 그래프 정렬을 변이 콜링까지 이어붙이기

PanVariants(bioRxiv, 2026)는 그래프 정렬 결과를 SNV·indel·SV·CNV·STR 콜링까지 이어붙인 프레임워크입니다. 이 프리프린트의 비교에서는 참조 패널의 구성과 크기, 시퀀싱 플랫폼, 변이 유형에 따라 재현율과 F1이 달라졌습니다. 예를 들어 논문이 보고한 한 비교에서는 참조 패널을 88개에서 214개 haplotype으로 확장했을 때 recall이 62.46%에서 68.28%, F1이 72.57%에서 75.23%로 상승했습니다. 따라서 특정 수치를 모든 데이터에 일반화하기보다, 선형 참조가 충분히 대표하지 못하는 변이를 그래프가 얼마나 회수하는지를 동일한 truth set과 평가 구간에서 확인해야 합니다. 아직 동료평가 전 프리프린트이므로 수치와 결론은 최종 출판본에서 다시 확인해야 합니다.

$\text{선형 파이프라인 변이 회수율} \;\lt\; \text{그래프 파이프라인 변이 회수율} \quad \text{(고다형성 영역일수록 격차 확대)}$

실습: 소규모 변이 그래프 구성과 그래프 정렬 (Colab CPU)

전체 인간 pangenome 그래프는 로컬에서 다루기 벅찬 크기이므로, 하나의 유전자 구간만 잘라 그래프의 동작 원리를 직접 확인합니다.

# Colab 또는 로컬 리눅스 셸에서 실행 (vg toolkit 정적 바이너리 사용)

# 1. vg 바이너리 받기 (릴리스 페이지의 정적 빌드)
wget https://github.com/vgteam/vg/releases/latest/download/vg -O vg
chmod +x vg

# 2. 예시 참조 서열(작은 FASTA)과 변이(VCF)로 그래프 구성
#    ref.fa: 임의의 1kb 구간, sites.vcf: 그 구간 안 SNP·indel 몇 개
./vg construct -r ref.fa -v sites.vcf.gz > graph.vg

# 3. 그래프 인덱싱 — 미니마이저 시딩용 인덱스와 GBWT 하플로타입 인덱스
./vg autoindex --workflow giraffe -r ref.fa -v sites.vcf.gz -p graph_index

# 4. 시뮬레이션 리드로 그래프 정렬
./vg giraffe -Z graph_index.giraffe.gbz -m graph_index.min -d graph_index.dist \
  -f reads.fastq > aligned.gam

# 5. 정렬 통계 확인 — 참조 위치가 아니라 그래프 경로 기준 정렬률
./vg stats -a aligned.gam

vg stats로 나오는 정렬률을 같은 리드를 BWA-MEM2로 선형 정렬했을 때와 비교해보면, 변이가 몰려 있는 구간일수록 그래프 정렬의 정렬률이 더 높게 나오는 경향을 직접 확인할 수 있습니다. 실제 임상·연구 파이프라인에서는 이 단계 뒤에 그래프 좌표를 다시 GRCh38 선형 좌표로 투영(surject)하는 과정이 이어지는데, 이 좌표 변환 자체가 pangenome 도입의 실무적 걸림돌 중 하나입니다.

CS 매핑

• bidirected sequence graph: variation graph는 노드가 서열 조각, 간선이 조각 끝 사이의 인접 관계를 나타내며 역위나 반복 때문에 순환을 포함할 수 있습니다. 따라서 본질적으로 DAG라고 볼 수 없고, 정렬은 이 그래프에서 리드와 잘 맞는 경로를 찾는 문제로 재정의됩니다.
• succinct 자료구조: GBWT는 그래프에 담긴 수십만 개 하플로타입 경로를 원본 대비 훨씬 작은 공간에 압축해 저장하면서도 빠른 질의를 지원하는 succinct 자료구조로, M14～M15의 BWT·FM-index를 그래프·다중 서열 맥락으로 확장한 것입니다.
• 그래프 탐색 vs 선형 탐색: 선형 정렬의 seed-and-extend가 좌표 하나를 기준으로 확장한다면, 그래프 정렬은 여러 경로 후보를 동시에 탐색하는 탐색 공간 확장 문제와 같은 구조를 가집니다.

자주 만나는 결함

• "pangenome이 GRCh38을 완전히 대체한다"는 오해: 대부분의 실무 파이프라인은 여전히 GRCh38(또는 T2T-CHM13) 좌표계를 결과 보고·데이터베이스 연동의 기준으로 삼습니다. 그래프 정렬 결과를 이 선형 좌표계로 투영하는 과정이 필요하며, 그래프만으로 기존 좌표계 기반 도구·주석을 전부 대체하지는 못합니다.
• 그래프가 커질수록 복잡도가 선형이 아니라는 점 간과: 구조 변이가 밀집한 구간은 그래프가 뒤엉킨 서브그래프(tangle)를 형성해 인덱싱·정렬 비용이 예상보다 크게 늘어날 수 있습니다. 전체 게놈을 무작정 그래프화하기 전에, 어떤 구간을 그래프로 표현할지 범위를 정하는 편이 안전합니다.

더 깊게 파고 싶다면

본문은 BPD 연구진이 직접 재구성한 서술입니다. 원 논문과 공식 자료로 심화해봅시다.

• Human Pangenome Reference Consortium (Liao et al., 2023, Nature): 47명 94 하플로타입 초안 범유전체 원 논문.
• vg toolkit 공식 GitHub (vgteam/vg): 그래프 구성·인덱싱·정렬 전체 명령어 레퍼런스.
• GBWT 원 논문 (Sirén et al.): 그래프 위 하플로타입 압축 인덱싱 설계.
• PanVariants (bioRxiv, 2026): 그래프 정렬–변이 콜링 통합 파이프라인. 정식 출판 시 인용 정보 갱신 필요.

F2 소주제의 "미래 편" 묶음(F32～F35) 중 첫 편으로, 참조 유전체 자체를 바꾸는 접근을 다뤘습니다. 다음 편 F33에서는 이 그래프·파운데이션 모델 인프라 위에서 움직이는 AI 에이전트로 시야를 넓힙니다.