4강 / 20강

VCF 여덟 열 — 파일을 실제로 읽는다

CHROM POS ID REF ALT QUAL FILTER INFO

VCF 는 탭으로 나뉜 표입니다. ## 로 시작하는 머리말이 길게 붙고, #CHROM 줄부터가 표의 열 이름입니다.

##fileformat=VCFv4.2

##reference=file://kit_locus.fa

##INFO=<ID=DP,Number=1,Type=Integer,Description="Raw read depth">

##FORMAT=<ID=AD,Number=R,Type=Integer,Description="Allelic depths">

#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT tumor

13 47258004 . A T 228.4 PASS DP=61;MQ=60 GT:AD:DP 1/1:10,51:61

앞 여덟 열은 어느 파일에나 있습니다. 아홉째부터가 시료마다 붙는 열입니다.

#	열	무엇	함정
1	CHROM	염색체 이름	`13` 인가 `chr13` 인가 (I 권 7강)
2	POS	자리. 1부터, 끝 포함	BED 는 0부터다. 한 칸 어긋난다
3	ID	알려진 변이 번호 (rs…)	개는 대부분 `.` — 목록이 빈약하다
4	REF	참조의 염기	참조를 바꾸면 이것도 바뀐다
5	ALT	다른 염기. 쉼표로 여럿	다중대립이면 뒤의 AD·GT 해석이 달라진다
6	QUAL	이 자리에 변이가 있다는 확신 (Phred)	변이의 세기가 아니다. 깊이가 깊으면 그냥 커진다
7	FILTER	`PASS` 또는 걸린 이유	`.` 는 '필터를 안 돌렸다' 이지 PASS 가 아니다
8	INFO	자리 전체에 대한 주석. `key=value;…`	도구마다 키가 다르다

VCF 를 다루는 코드는 결국 이 몇 줄입니다. pysam·cyvcf2 가 이 일을 빠르게 해 줄 뿐입니다.

Python 브라우저에서 실행

# VCF 한 줄을 손으로 파싱한다
line = "13\t47258004\t.\tA\tT\t228.4\tPASS\tDP=61;MQ=60;VDB=0.42\tGT:AD:DP\t1/1:10,51:61"
cols = line.split("\t")
names = ["CHROM", "POS", "ID", "REF", "ALT", "QUAL", "FILTER", "INFO", "FORMAT", "tumor"]

for n, v in zip(names, cols):
    print(f"  {n:<8}{v}")
print()

info = dict(kv.split("=", 1) if "=" in kv else (kv, True) for kv in cols[7].split(";"))
print("INFO 를 풀면:", info)
print()

fmt = cols[8].split(":")
val = cols[9].split(":")
sample = dict(zip(fmt, val))
print("시료 열을 풀면:", sample)
print()

ad = [int(x) for x in sample["AD"].split(",")]
dp = int(sample["DP"])
print(f"참조 리드 {ad[0]}   변이 리드 {ad[1]}   깊이 {dp}")
print(f"VAF = {ad[1]}/{dp} = {ad[1]/dp:.3f}")

FILTER 열의 . 는 PASS 가 아니다

가장 자주 나는 조용한 사고입니다.

값	뜻
`PASS`	필터를 돌렸고 통과했다
`.`	필터를 아직 안 돌렸다. 통과했다는 뜻이 아니다
`LowQual;strand_bias`	이 이유들로 걸렸다

bcftools call 의 날 출력은 FILTER 가 전부 . 입니다. 그것을 grep PASS 로 거르면 0줄이 나오고, grep -v FAIL 로 거르면 전부가 통과합니다.

워크숍 노트북은 필터를 손으로 겁니다 — 11강에서 그 코드를 봅니다.

QUAL 문턱을 하나로 정하면 얕은 자리가 통째로 사라집니다. 커버리지가 고르지 않을 때 특히 아픕니다.

Python 브라우저에서 실행

# QUAL 이 무엇을 재는가 — 깊이가 깊으면 그냥 커진다
from scipy.stats import binomtest
import math

err = 0.002
print(f"{'상황':<24}{'AD':>10}{'VAF':>8}{'QUAL':>9}")
for name, alt, dp in [("얕은 이형접합", 5, 10), ("깊은 이형접합", 30, 60),
                      ("아주 깊은 이형접합", 150, 300),
                      ("얕은 subclone", 2, 10), ("깊은 subclone", 12, 300)]:
    p = binomtest(alt, dp, err, alternative="greater").pvalue
    qual = min(999, -10 * math.log10(max(p, 1e-99)))
    print(f"{name:<24}{f'{dp-alt},{alt}':>10}{alt/dp:>8.3f}{qual:>9.1f}")
print()
print("맨 위 세 줄은 VAF 가 같은데 QUAL 이 다섯 배씩 벌어집니다.")
print("QUAL 은 '변이가 얼마나 센가' 가 아니라 '변이가 있다는 것을 얼마나 확신하는가' 입니다.")
print()
print("그래서 QUAL 로 줄 세우면 깊이 깊은 자리가 위로 올라옵니다.")
print("생물학적 중요도로 줄 세우려면 VAF 와 주석이 필요합니다 — 12·13강.")

개에서 ID 열이 대부분 `.` 인 것

사람 VCF 를 열면 ID 열에 rs334 같은 번호가 빽빽합니다. dbSNP 에 등록된 변이라는 뜻입니다.

개 VCF 는 대부분 . 입니다. 등록된 목록이 훨씬 작기 때문입니다.

🔍

알려진 변이 확인

"이건 흔한 다형성이다" 를 못 한다

🧪

BQSR

I 권 18강 — 학습 자료가 오염된다

➖

생식세포 걸러내기

gnomAD 같은 빈도표가 없어 혈액에만 기댄다

세 번째가 이 권에서 중요합니다. 사람이라면 혈액이 없어도 집단 빈도표로 "이건 사람 1%가 갖고 있는 다형성이다" 를 걸러 낼 수 있습니다. 개는 그것이 안 되므로 짝을 맞춰 읽는 것 말고는 방법이 없습니다.

3,000 호주달러의 절반을 병 없는 조직에 쓴 것이 개에서는 선택이 아니라 필수였습니다.

🤔 확인 문제

VCF 의 FILTER 열이 . 이면?

해설 bcftools call 의 날 출력은 FILTER 가 전부 . 입니다. 이것을 grep PASS 로 거르면 0줄이 나오고, 반대로 '걸린 것만 뺀다'로 거르면 전부 통과합니다. 두 사고가 정반대 방향으로 조용히 납니다.

🤔 확인 문제

QUAL 열이 재는 것은?

🤔 확인 문제

개 VCF 의 ID 열이 대부분 . 라는 사실이 이 권에서 갖는 가장 실질적인 결과는?

개에서 ID 열이 대부분 . 인 것

개에서 ID 열이 대부분 `.` 인 것