14강 / 20강

MAPQ — 확신을 숫자로

제6장 「그 자리를 얼마나 확신하는지도 함께 적어 준다」

소설의 이 문장이 MAPQ 의 정의입니다.

"거기에 떨어진 조각은 어디서 왔는지 알 수 없다. 그래서 정렬 프로그램은 자리 하나만 정해 주는 것이 아니라 그 자리를 얼마나 확신하는지도 함께 적어 준다. 확신이 낮은 조각은 나중에 변이를 찾을 때 무게를 낮추거나 버린다."

MAPQ = −10 · log₁₀ P(자리가 틀렸을 확률)

Phred 와 같은 자입니다. 다만 재는 것이 다릅니다 — 염기가 틀렸을 확률이 아니라 자리가 틀렸을 확률입니다.

MAPQ 는 표준화되어 있지 않습니다. 규격이 정의만 정해 두었고 계산은 도구가 알아서 합니다.

Python 브라우저에서 실행

# MAPQ 를 읽는 표
def p_wrong(q):
    return 10 ** (-q / 10)

rows = [
    (0,  "어디인지 모른다 — 똑같이 좋은 자리가 여럿"),
    (1,  "거의 모른다"),
    (3,  "후보가 둘쯤"),
    (10, "10 번에 한 번은 틀린다"),
    (20, "100 번에 한 번"),
    (30, "1,000 번에 한 번"),
    (60, "BWA-MEM 의 최대값 — 사실상 유일한 자리"),
]
print(f"{'MAPQ':>6}{'자리가 틀릴 확률':>18}   뜻")
for q, memo in rows:
    print(f"{q:>6}{p_wrong(q):>18.4f}   {memo}")
print()
print("주의: 도구마다 최대값이 다릅니다.")
print("  BWA-MEM   60")
print("  Bowtie2   42")
print("  minimap2  60")
print("그래서 'MAPQ >= 30' 같은 필터를 도구를 바꾸며 그대로 쓰면 뜻이 달라집니다.")

MAPQ 0 이 뜻하는 것 — 가장 자주 오해받는 값

흔한 오해

"이 리드는 붙지 않았다" · "품질이 나쁘다"

실제

"똑같이 좋은 자리가 여럿이라 하나를 고를 수 없다"

리드 자체는 멀쩡할 수 있습니다. 참조 유전체가 그 서열을 여러 군데에 갖고 있는 것뿐입니다. 정렬기는 그중 하나를 임의로 골라 적고 MAPQ 0 을 붙입니다.

그래서 MAPQ 0 리드를 그냥 두면 가짜 변이가 생깁니다 — 다른 자리에서 온 리드가 여기 있는 것처럼 보이니까요. 거의 모든 변이 부르기가 -q 10 이나 -q 20 을 겁니다.

실제 BWA-MEM 공식은 이보다 복잡하지만, 1등과 2등의 차이라는 뼈대는 같습니다.

Python 브라우저에서 실행

# MAPQ 가 어떻게 계산되나 — BWA-MEM 의 대강
def mapq_like(best, second, is_unique=False):
    "1등과 2등 정렬 점수의 차이가 확신이 된다 (단순화한 모형)"
    if is_unique:
        return 60
    if second is None:
        return 60
    gap = best - second
    return max(0, min(60, int(gap * 2.5)))

cases = [
    ("반복 없는 구간",     150, None, True),
    ("작은 차이가 있는 중복", 150, 145, False),
    ("거의 같은 중복",     150, 148, False),
    ("완전히 같은 중복",   150, 150, False),
    ("SINE 반복 안",       120, 119, False),
]
print(f"{'상황':<22}{'1등':>6}{'2등':>6}{'MAPQ':>7}   해석")
for name, b, s, uq in cases:
    q = mapq_like(b, s, uq)
    verdict = "믿는다" if q >= 30 else ("애매" if q > 0 else "자리를 모른다")
    print(f"{name:<22}{b:>6}{str(s):>6}{q:>7}   {verdict}")
print()
print("요점은 절대 점수가 아니라 1등과 2등의 '차이' 입니다.")
print("아주 잘 붙었어도 똑같이 잘 붙는 자리가 또 있으면 MAPQ 는 0 입니다.")

개 유전체에서 이것이 왜 더 아픈가

🕳️

메우지 못한 구멍

CanFam3.1 에 23,000개 이상. 거기 살던 리드가 비슷한 다른 자리로 간다

📋

덜 정확한 주석

반복 구간 목록(RepeatMasker 트랙)이 사람 것만큼 잘 정리돼 있지 않다

🐕

품종 거리

로지는 믹스견, 참조는 복서. 차이가 큰 자리에서 정렬 점수가 함께 낮아진다

세 번째가 특히 고약합니다. 참조와 먼 자리는 1등 점수도 낮아지는데, 그러면 2등과의 차이도 함께 줄어 MAPQ 가 떨어집니다. "내 개가 참조와 다르다"가 "여기는 못 믿겠다"로 번역되는 것입니다.

이 곡선의 어깨가 어디인지는 데이터마다 다릅니다. 그래서 문턱을 정하기 전에 한 번 그려 봅니다.

Python 브라우저에서 실행

# 필터를 옮기면 무엇이 남고 무엇이 사라지나
import numpy as np
rng = np.random.default_rng(7)

n = 20000
# 대부분은 유일하게 붙고(60), 일부는 반복 구간에 있다
mapq = np.concatenate([
    np.full(int(n * .78), 60),
    rng.integers(25, 60, int(n * .10)),
    rng.integers(1, 25,  int(n * .07)),
    np.zeros(int(n * .05), dtype=int),
])
rng.shuffle(mapq)

print(f"{'문턱':>6}{'남는 리드':>12}{'비율':>9}")
for thr in (0, 1, 10, 20, 30, 60):
    keep = (mapq >= thr).sum()
    print(f"{thr:>6}{keep:>12,}{keep/len(mapq):>9.1%}")
print()
print("0 → 1 에서 5% 가 사라집니다. 이것이 '자리를 모르는' 리드입니다.")
print("20 → 30 에서 또 몇 % 가 사라집니다. 여기부터는 진짜 신호도 함께 잘립니다.")

🤔 확인 문제

MAPQ 0 이 뜻하는 것은?

🤔 확인 문제

MAPQ 를 결정하는 핵심 요소는?

🤔 확인 문제

개 데이터에서 MAPQ 필터를 사람과 똑같이 세게(-q 30) 걸면 생기는 문제는?

해설 로지는 믹스견이고 참조는 복서 타샤입니다. 참조와 다른 자리는 1등 점수가 낮아지고 2등과의 차이도 줄어 MAPQ 가 떨어집니다. '내 개가 참조와 다르다'가 '여기는 못 믿겠다'로 번역되는 것입니다. 실무는 -q 10~20 을 쓰고 낮은 자리에 꼬리표를 답니다.