BioPlayground

🧬
목록으로

FASTQ 포맷과 품질 점수: Phred33/64 인코딩의 실무 감각

NGS 실무의 가장 첫 파일인 FASTQ. 네 줄이 어떻게 얽혀 있는지, Phred 점수가 왜 로그 스케일인지, Illumina가 왜 인코딩을 두 번 갈아엎었는지, 파이썬 20줄 스트리밍 파서로 손에 익혀봅시다.

입문
|
15
|
검증 완료 (2026-07-22)
NGSsequencing qualityFASTQ format
진행률0/52 (0%)

"FASTQ to Paper"를 시작하며

여기서부터가 "FASTQ to Paper" 시리즈의 실질적 출발점입니다. 앞의 Micro 티어(M01~M30)에서 우리는 알고리즘의 알파벳 — 정렬·인덱싱·HMM·계통수 — 을 손으로 유도했습니다. 이제부터 25편에 걸쳐, 그 알고리즘이 어떤 실무 파일 포맷을 타고, 어떤 도구를 지나, 어떤 논문 그림으로 마무리되는지를 따라갑니다.

그리고 이 여정의 첫 파일이 바로 FASTQ입니다. 시퀀서에서 나온 원시 리드 하나하나가 이 포맷에 담겨 있고, 이 파일이 잘못되면 뒤의 모든 파이프라인이 뒤틀립니다. 그러니 서두르지 맙시다.

네 줄이 한 리드입니다

FASTQ 파일 한 리드는 정확히 네 줄입니다.

text
@SRR1039508.1 HWI-ST177:290:C0TECACXX:1:1101:1225:2130/1
CATTGCTGATACCAATTAGGTCGTGAGGACTGGGCAAGTGTAAAGCT
+
HHHHHHHHHHHHHHHHHHHHHHHHHHHHHGGFHHHHHHHGGGGHHFB
  • 1행 (헤더): @로 시작하는 리드 식별자. 뒤의 문자열은 시퀀서·플로우셀·타일·좌표를 인코딩합니다.
  • 2행 (서열): A/C/G/T/N 문자열. 이 리드의 실제 염기 서열.
  • 3행 (구분자): +. 종종 헤더가 반복되기도 하지만 실무에선 그냥 + 한 글자입니다.
  • 4행 (품질 문자열): 서열과 정확히 같은 길이의 ASCII 문자열. 각 문자 하나가 해당 위치 염기의 품질 점수입니다.

이 네 줄이 파일 안에 수백만~수억 개 반복되어 있습니다. 30배 커버리지 인간 게놈이면 리드 1억 개, 파일 크기는 gzip 압축으로도 30~50GB급입니다.

Phred 점수: 로그로 표현한 오류 확률

품질 문자열이 왜 그냥 HHHHHH... 같은 알파벳으로 되어 있을까요. 그 문자 하나가 실은 Phred 점수라는 정수를 ASCII로 인코딩한 것입니다.

Phred 점수 Q는 이렇게 정의됩니다.

Q=10log10(Perror)Q = -10 \cdot \log_{10}(P_{\text{error}})
  • Q10 → 오류 확률 0.1 (10 리드 중 1개 틀림)
  • Q20 → 오류 확률 0.01
  • Q30 → 오류 확률 0.001
  • Q40 → 오류 확률 0.0001

실무 감각으로 외워둘 값이 있습니다. **Q30이 "게놈 재분석 가능한 하한선"**입니다. Illumina NovaSeq가 뽑는 리드 대부분이 Q30 이상이고, PacBio HiFi도 Q20~Q40 사이에 걸칩니다. Q10 이하가 다수인 리드가 있다면 그건 어댑터 잔재 또는 시약 소진 신호로 봐야 합니다.

로그를 왜 씌웠는가는 계산의 편의 때문입니다. 확률 0.0001을 파일에 그대로 담으면 4자리인데, Q40이면 정수 두 자리로 끝납니다. 그리고 여러 위치 확률을 곱해야 할 때, 로그 스케일에서는 그냥 더하면 됩니다. 이건 CS의 로그 확률 트릭 그대로입니다.

Phred33 vs Phred64 — 두 인코딩의 역사

이제 진짜 함정입니다. 정수 Q를 ASCII로 어떻게 매핑하느냐가 두 가지 규격으로 갈립니다.

인코딩시작 문자ASCII 시작값범위시대
Phred+33 (Sanger, Illumina 1.8+)!33Q0 ~ Q41 (! ~ J)2011~ 현재 표준
Phred+64 (Illumina 1.3 ~ 1.7)@64Q0 ~ Q41 (@ ~ h)2004~2011 구형

즉 파일에서 문자 H를 봤다면,

  • Phred33이면 Q = ord('H') − 33 = 72 − 33 = Q39 (거의 완벽)
  • Phred64이면 Q = ord('H') − 64 = 72 − 64 = Q8 (거의 오류)

같은 문자가 두 세계에서 완전히 다른 의미를 갖는다는 뜻입니다. 오래된 SRA 데이터를 내려받아 그대로 파이프라인에 넣었다가 "품질이 왜 이 모양이지?"로 이틀을 태우는 사고가 여기서 옵니다.

자동 감지의 원칙

파일이 어느 인코딩인지 어떻게 알까요. 실무 규칙은 이렇습니다.

  1. 첫 몇만 리드의 품질 문자 최소값을 봅니다.
  2. 최소값이 ASCII 33~58 사이라면 → Phred33 (58 이하 문자는 Phred64에 없기 때문).
  3. 최소값이 ASCII 64 이상만 나타난다면 → Phred64 의심.
  4. 애매하면 데이터셋의 발표 연도로 결정 (2012년 이후 = 거의 확실히 Phred33).

FastQC가 처음 실행될 때 콘솔에 Encoding: Sanger / Illumina 1.9이라고 뜨는 게 정확히 이 판정입니다.

파이썬 스트리밍 파서 20줄

파일이 30GB짜리라 한 번에 메모리에 못 올립니다. 스트리밍으로 네 줄씩 뜯는 게 실무 원칙입니다.

python
import gzip
def iter_fastq(path):
opener = gzip.open if path.endswith(".gz") else open
with opener(path, "rt") as f:
while True:
header = f.readline()
if not header:
return
seq = f.readline().rstrip()
plus = f.readline()
qual = f.readline().rstrip()
yield header.rstrip(), seq, qual
def phred33_to_q(qual_str):
return [ord(c) - 33 for c in qual_str]
# 사용 예: 첫 100 리드의 평균 품질
for i, (hdr, seq, qual) in enumerate(iter_fastq("reads.fastq.gz")):
if i >= 100:
break
qs = phred33_to_q(qual)
print(f"read {i}: len={len(seq)}, meanQ={sum(qs)/len(qs):.1f}")

이게 전부입니다. gzip 여부 자동 판단, 스트리밍, ASCII 33 감산. 파일이 3TB짜리여도 이 루프는 상수 메모리로 돕니다.

서열 길이 ≠ 품질 문자열 길이 검증

FASTQ 파서를 처음 짤 때 실무에서 반드시 넣는 방어 코드는 이겁니다.

python
for hdr, seq, qual in iter_fastq(path):
if len(seq) != len(qual):
raise ValueError(f"length mismatch: {hdr}")

시퀀서 오류·전송 손실·부분 gzip 잘림 같은 사건들이 이 assertion에서 잡힙니다. 파이프라인 뒷단(BWA·GATK)이 조용히 죽어가는 것보다 여기서 크게 터지는 편이 훨씬 낫습니다.

실무 감각 — 세 가지 지표

FastQC 리포트를 열면 20개 이상의 항목이 뜹니다. 처음엔 다 볼 필요 없이 세 가지만 봅니다.

  1. Per base sequence quality: 좌→우 진행하면서 Q30 상자가 얼마나 오래 초록색을 유지하는지. 리드 뒷부분이 Q20 이하로 무너지면 3′ 트리밍이 필요합니다.
  2. Overrepresented sequences: 어댑터 잔재 검출. 특정 서열이 리드의 5% 이상 반복되면 어댑터/PCR 프라이머 오염 신호입니다.
  3. Sequence duplication levels: PCR 증폭 편향 신호. RNA-seq에선 자연 중복이 있지만, WGS에서 30% 이상 duplication은 라이브러리 재제작 결재감입니다.

다음 편(S02)에서 FastQC와 자동 트리밍 도구(fastp · Trim Galore)로 이 셋을 실제로 해결해봅니다.

Colab에서 실전 파일 열어보기

무료 환경에서 실제 FASTQ 감을 잡는 가장 빠른 길은 SRA의 소형 샘플입니다. Colab에서 아래 셀 두 개만 실행해보세요.

bash
!pip install pyfastx -q
!wget -q https://sra-pub-src-1.s3.amazonaws.com/SRR1039508/SRR1039508_1.fastq.gz -O r1.fq.gz
python
import pyfastx
fq = pyfastx.Fastq("r1.fq.gz")
print(f"총 리드: {len(fq):,}, 평균 리드 길이: {fq.avglen:.1f}bp")
print(f"GC%: {fq.gc_content:.1f}, 인코딩: {fq.phred}")
# 첫 3 리드 품질 요약
for name, seq, qual in fq[:3]:
qs = [ord(c) - fq.phred for c in qual]
print(name, "meanQ=", sum(qs) / len(qs))

fq.phred가 33이나 64로 뜨는 순간이 인코딩 자동 감지가 성공한 지점입니다.

각주 (실습 환경): Colab 무료 T4 티어면 충분합니다. Colab 한국 결제 시 우편번호에 미국 형식(예: 94043)을 넣으면 결제 화면을 통과할 수 있습니다. SRA 다운로드가 느리면 EBI ENA의 미러(https://www.ebi.ac.uk/ena/browser/)를 사용해도 됩니다.

CS 매핑

FASTQ는 사실상 세 가지 CS 개념의 합주입니다.

  • ASCII 오프셋 인코딩: 정수를 문자로 담기 위해 시작점을 옮긴다. Phred33/64는 오프셋 33·64의 관례일 뿐입니다.
  • 로그 확률(log-probability): 확률을 곱할 일이 많으니 로그 스케일로 저장. Q = −10log₁₀P.
  • 스트리밍 파서: 파일이 메모리를 초과하므로 라인 단위 이터레이터로 처리. Unix pipe 철학과 동일합니다.

DryBench의 "ASCII 인코딩과 텍스트 파싱" 편(ascii-encoding-and-text-parsing)이 이 개념들의 순수 CS 판입니다. FASTQ는 그 개념이 바이오 실무로 육화된 첫 사례라고 봐도 됩니다.

마무리 — 다음 편으로

지금 우리는 파일을 열 수 있고, 리드마다 품질 벡터를 뽑을 수 있으며, 인코딩을 자동 감지할 수 있습니다. 다음 편 S02에서는 이 리드에서 어댑터를 잘라내고, 3′ 저품질 꼬리를 다듬어 파이프라인 입구를 정돈합니다. 그 뒤 S03에서 BWA-MEM으로 참조 게놈에 실제로 정렬합니다.

더 깊게 파고 싶다면

본문은 BPD 자체 재구성입니다. 원리를 확실히 잡고 싶다면 아래 자료를 큐레이션합니다.

  • Illumina Knowledge Base — Understanding Illumina Quality Scores (공식 백서 링크). Q 점수 캘리브레이션 실측 그래프가 인상적입니다.
  • NCBI SRA Handbook의 FASTQ 섹션. 인코딩 변천의 공식 연대기.
  • Ewing & Green (1998), Base-calling of automated sequencer traces using phred. II. Genome Research 8:186. Phred 알고리즘 원 논문 — 오늘날 우리가 쓰는 Q 점수의 출발.
  • Broad Institute BroadE 유튜브 — Illumina Chemistry 15분 강의. 왜 3′ 품질이 떨어지는지 광학·화학적 원인이 시각화됩니다.

파일을 열어보고 손으로 파서를 짜본 사람만이 다음 편의 트리밍 결정에 자신감을 가질 수 있습니다. Colab을 켜고 SRR 하나 내려받아 봅시다.