11. 맞춰 놓기 — 150글자를 24억 칸에 · 손박사의 애완견 비만세포암 치료제 첫걸음 I

소설이 쓴 비유가 정확합니다.

"책 한 권을 잘게 찢어 종잇조각 수억 개를 만들었다고 하자. 조각 하나에는 글자 백오십 개가 적혀 있다. 이제 그 조각들이 원래 책의 어느 쪽 어느 줄에서 나왔는지 알아내야 한다. 다행히 그 책의 다른 판본이 한 권 있다. 판본이 참조 게놈이다."

그래서 퍼즐보다는 교정에 가깝습니다. 조각이 들어갈 자리는 이미 정해져 있고, 찾기만 하면 됩니다.

이론상 20글자면 충분합니다. 실제로는 반복 서열 때문에 150글자로도 모자란 자리가 있습니다.

Python 브라우저에서 실행

# 왜 150 글자면 대체로 자리가 하나로 정해지나 — 확률로
G = 2_400_000_000

print("길이 k 짜리 서열이 24억 글자 유전체에 우연히 또 나타날 기대 횟수")
print()
for k in (10, 15, 20, 25, 30, 50, 100, 150):
    expected = G / (4 ** k)
    if expected >= 1:
        print(f"  k = {k:>4}   {expected:>18,.0f} 번   ← 여러 군데")
    else:
        print(f"  k = {k:>4}   {expected:>18.2e} 번")
print()
print("k = 20 을 넘으면 우연한 일치가 사실상 사라집니다.")
print("그런데 실제 유전체는 무작위가 아닙니다 — 다음 셀.")

씨앗 길이를 바꿔 보세요. 짧으면 후보가 폭발하고, 길면 변이가 있는 리드를 놓칩니다.

불러오는 중…

정렬기가 실제로 하는 일 — 씨앗과 확장

🌱

1. 씨앗 찾기

리드에서 딱 맞는 짧은 토막을 색인으로 즉시 찾는다

↔️

2. 확장

그 자리에서 양쪽으로 늘려 가며 점수를 매긴다 (Smith-Waterman)

🏅

3. 고르기

가장 높은 점수 하나. 2등과의 차이가 MAPQ 가 된다

1단계가 빠른 이유는 FM-index (BWT) 라는 자료구조 덕분입니다. 참조 유전체를 한 번 뒤집어 정리해 두면, 임의의 짧은 문자열이 어디에 있는지를 유전체 길이와 무관하게 찾을 수 있습니다.

그래서 bwa index 를 한 번 돌려 두면 그다음부터는 빠릅니다. 코드랩 3단계의 bwa index kit_locus.fa 가 그것입니다.

정렬기 매개변수는 대부분 이런 맞바꿈입니다. 어느 쪽 실수를 감수할지 고르는 것뿐입니다.

Python 브라우저에서 실행

# 씨앗 길이의 딜레마 — 짧으면 폭발하고 길면 놓친다
G = 2_400_000_000
read_len = 150

print(f"{'씨앗 k':>7}{'우연 후보 수':>16}{'변이 1개 있는 리드에서':>26}")
print(f"{'':>7}{'(자리당)':>16}{'k-mer 하나라도 온전할 확률':>26}")
print()
for k in (11, 15, 19, 25, 31, 51):
    noise = G / (4 ** k)
    # 리드 안에 변이가 하나 있을 때, 그것을 비껴가는 온전한 k-mer 가 존재할 확률
    clean = 1.0 if read_len - k >= k else 0.0
    n_windows = read_len - k + 1
    n_clean = max(0, n_windows - k)          # 변이를 덮지 않는 창의 수
    print(f"{k:>7}{noise:>16.2e}{n_clean:>20} 개")
print()
print("k 를 키우면 잡음은 줄지만, 변이·오류를 비껴가는 창도 줄어듭니다.")
print("BWA-MEM 의 기본 씨앗은 19 입니다. 그 사이 어디쯤이 좋다는 경험값입니다.")

20 kb 짜리 교실 참조의 한계 — 미리 적어 둡니다

워크숍 노트북은 KIT 유전자좌 20 kb 만 참조로 씁니다. 그러면 정렬이 아주 쉬워집니다 — 후보 자리가 24억 칸이 아니라 2만 칸뿐이니까요.

그래서 얻는 것

교실에서 몇 초 만에 돈다. 논리를 눈으로 따라갈 수 있다

그래서 잃는 것

반복 구간이 거의 없다. MAPQ 가 전부 높게 나온다. 실제 WGS 의 어려움이 사라진다

그래서 이 강의는 MAPQ 를 따로 14강에서 다룹니다. 노트북만 돌려서는 배울 수 없는 것이기 때문입니다.

🤔 확인 문제

정렬기가 리드를 참조에 붙이는 순서로 옳은 것은?

해설 씨앗(seed) → 확장(extend) → 선택입니다. 씨앗 찾기가 빠른 것은 FM-index(BWT) 덕분으로, 유전체 길이와 거의 무관하게 짧은 문자열의 위치를 찾습니다. 그래서 bwa index 를 한 번 돌려 두는 것입니다.

🤔 확인 문제

씨앗 길이 k 를 크게 잡으면 무슨 일이 생깁니까?

🤔 확인 문제

워크숍 노트북이 20 kb 짜리 참조를 쓰기 때문에 배울 수 없게 되는 것은?