"FASTQ to Paper"를 시작하며
여기서부터가 "FASTQ to Paper" 시리즈의 실질적 출발점입니다. 앞의 Micro 티어(M01~M30)에서 우리는 알고리즘의 알파벳 — 정렬·인덱싱·HMM·계통수 — 을 손으로 유도했습니다. 이제부터 25편에 걸쳐, 그 알고리즘이 어떤 실무 파일 포맷을 타고, 어떤 도구를 지나, 어떤 논문 그림으로 마무리되는지를 따라갑니다.
그리고 이 여정의 첫 파일이 바로 FASTQ입니다. 시퀀서에서 나온 원시 리드 하나하나가 이 포맷에 담겨 있고, 이 파일이 잘못되면 뒤의 모든 파이프라인이 뒤틀립니다. 그러니 서두르지 맙시다.
네 줄이 한 리드입니다
FASTQ 파일 한 리드는 정확히 네 줄입니다.
@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는 이렇게 정의됩니다.
- 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+) | ! | 33 | Q0 ~ Q41 (! ~ J) | 2011~ 현재 표준 |
| Phred+64 (Illumina 1.3 ~ 1.7) | @ | 64 | Q0 ~ Q41 (@ ~ h) | 2004~2011 구형 |
즉 파일에서 문자 H를 봤다면,
- Phred33이면 Q = ord('H') − 33 = 72 − 33 = Q39 (거의 완벽)
- Phred64이면 Q = ord('H') − 64 = 72 − 64 = Q8 (거의 오류)
같은 문자가 두 세계에서 완전히 다른 의미를 갖는다는 뜻입니다. 오래된 SRA 데이터를 내려받아 그대로 파이프라인에 넣었다가 "품질이 왜 이 모양이지?"로 이틀을 태우는 사고가 여기서 옵니다.
자동 감지의 원칙
파일이 어느 인코딩인지 어떻게 알까요. 실무 규칙은 이렇습니다.
- 첫 몇만 리드의 품질 문자 최소값을 봅니다.
- 최소값이 ASCII 33~58 사이라면 → Phred33 (58 이하 문자는 Phred64에 없기 때문).
- 최소값이 ASCII 64 이상만 나타난다면 → Phred64 의심.
- 애매하면 데이터셋의 발표 연도로 결정 (2012년 이후 = 거의 확실히 Phred33).
FastQC가 처음 실행될 때 콘솔에 Encoding: Sanger / Illumina 1.9이라고 뜨는 게 정확히 이 판정입니다.
파이썬 스트리밍 파서 20줄
파일이 30GB짜리라 한 번에 메모리에 못 올립니다. 스트리밍으로 네 줄씩 뜯는 게 실무 원칙입니다.
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 파서를 처음 짤 때 실무에서 반드시 넣는 방어 코드는 이겁니다.
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개 이상의 항목이 뜹니다. 처음엔 다 볼 필요 없이 세 가지만 봅니다.
- Per base sequence quality: 좌→우 진행하면서 Q30 상자가 얼마나 오래 초록색을 유지하는지. 리드 뒷부분이 Q20 이하로 무너지면 3′ 트리밍이 필요합니다.
- Overrepresented sequences: 어댑터 잔재 검출. 특정 서열이 리드의 5% 이상 반복되면 어댑터/PCR 프라이머 오염 신호입니다.
- Sequence duplication levels: PCR 증폭 편향 신호. RNA-seq에선 자연 중복이 있지만, WGS에서 30% 이상 duplication은 라이브러리 재제작 결재감입니다.
다음 편(S02)에서 FastQC와 자동 트리밍 도구(fastp · Trim Galore)로 이 셋을 실제로 해결해봅니다.
Colab에서 실전 파일 열어보기
무료 환경에서 실제 FASTQ 감을 잡는 가장 빠른 길은 SRA의 소형 샘플입니다. Colab에서 아래 셀 두 개만 실행해보세요.
!pip install pyfastx -q!wget -q https://sra-pub-src-1.s3.amazonaws.com/SRR1039508/SRR1039508_1.fastq.gz -O r1.fq.gzimport pyfastxfq = 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 하나 내려받아 봅시다.