BioPlayground

🧬
목록으로

samtools 실무: BAM/CRAM 흐름 · view · sort · index · CRAM 압축의 감각

samtools의 명령 하나하나가 무엇을 하는지, BAM에서 CRAM으로 넘어갈 때 저장 공간이 왜 절반 이하로 줄어드는지, 실무 파이프라인에서 어떤 순서로 부르는지 정리합니다.

입문
|
18분
|
검증 완료 (2026-07-22)
BAM manipulationCRAM compressionsamtools
진행률0/52 (0%)

samtools를 왜 익숙해져야 하나요

지난 두 편에서 BWA-MEM과 STAR/HISAT2로 정렬을 마쳤습니다. 이제 그 결과 BAM 파일을 다뤄야 합니다. GATK, IGV, 커스텀 파이썬 스크립트 — 이 모든 뒷단이 samtools의 명령들을 전제로 씁니다. samtools는 그 자체가 CLI 도구이자 htslib이라는 C 라이브러리의 얼굴입니다.

한 시간만 각 잡고 익히면 그 뒤 5년의 파이프라인 생활이 편해집니다. 서두르지 맙시다.

SAM · BAM · CRAM — 셋의 관계

세 포맷은 같은 정보를 다른 방식으로 저장합니다.

포맷인코딩크기용도
SAM텍스트 (탭 구분)원본, 매우 큼사람이 볼 때·디버깅·파이프 스트리밍
BAM이진 압축 (BGZF)SAM의 25～30%표준 저장·인덱스 지원·거의 모든 도구가 소비
CRAM레퍼런스 기반 압축BAM의 30～60%장기 저장·클라우드 아카이빙

전환은 자유롭습니다. samtools 하나로 세 방향 어디로든 갈 수 있습니다.

첫 명령 다섯 개면 충분합니다

실무에서 매일 부르는 명령은 다섯 개입니다.

1. samtools view — 파일을 열어봅시다

bash
samtools view SAMPLE01.bam | head -3
samtools view -h SAMPLE01.bam | head -20   # 헤더 포함
samtools view -c SAMPLE01.bam              # 총 리드 수
samtools view -c -f 4 SAMPLE01.bam          # unmapped 리드 수
samtools view -c -F 256 SAMPLE01.bam        # secondary 제외 리드 수
  • -h: 헤더도 출력.
  • -c: 개수만 세기.
  • -f N: 이 플래그가 켜진 리드만 (예: 4 = unmapped).
  • -F N: 이 플래그가 켜진 리드는 제외 (예: 256 = secondary).

2. samtools sort — 좌표 정렬 필수 관문

bash
samtools sort -@ 4 -o SAMPLE01.sorted.bam SAMPLE01.bam

BWA/STAR가 뱉는 SAM/BAM은 리드 순서 그대로입니다. GATK, IGV, bcftools — 대부분이 좌표 정렬(coordinate-sorted) BAM을 요구합니다. sort는 그 정렬을 수행합니다.

  • -@ 4: 스레드 4개.
  • -o out.bam: 출력.
  • 이름순 정렬이 필요하면 -n (특정 상황: 페어 재통합, HTSeq).

3. samtools index — 랜덤 액세스의 열쇠

bash
samtools index SAMPLE01.sorted.bam

.bai 파일이 생기고, 이후 IGV/bcftools/samtools view가 특정 구간 쿼리를 즉시 할 수 있습니다.

bash
samtools view SAMPLE01.sorted.bam chr17:41196312-41277500   # BRCA1 구간만

인덱스가 없으면 이 명령은 파일 전체를 스캔해서 몇 분이 걸립니다. 있으면 밀리초에 끝납니다. Snakemake 워크플로우에서 sort 뒤에 index를 반드시 이어붙이는 게 관습입니다.

4. samtools flagstat — 파이프라인 첫 대시보드

bash
samtools flagstat SAMPLE01.sorted.bam

출력 예시(요약):

text
120,000,000 + 0 in total
115,200,000 + 0 mapped (96.00% : N/A)
119,800,000 + 0 paired in sequencing
114,900,000 + 0 properly paired (95.91% : N/A)
2,400,000 + 0 duplicates

이 다섯 줄만 매일 봐도 파이프라인 건강 진단이 됩니다.

  • mapped ratio: 90% 아래면 참조 게놈 불일치, 오염 의심.
  • properly paired ratio: 85% 아래면 라이브러리 인서트 크기 문제 또는 어댑터 잔재.
  • duplicates: 20% 이상이면 라이브러리 복잡도 저하.

5. samtools idxstats — 염색체별 리드 분포

bash
samtools idxstats SAMPLE01.sorted.bam

각 시퀀스(염색체·컨티그) 별로 리드가 몇 개 붙었는지 알려줍니다. 이 결과로 성별 판정(chrY 리드 비율)이나 미토콘드리아 이상치를 즉시 확인합니다.

CRAM — 왜 실무가 넘어가고 있나요

BAM은 서열을 리드마다 다 저장합니다. 그런데 리드 대부분이 참조 게놈과 거의 같습니다. 참조와 다른 부분만 저장하면 훨씬 줄어들지 않을까? 이 아이디어를 구현한 게 CRAM입니다.

  • CRAM은 참조를 기준으로 각 리드가 어떤 차이가 있는지만 이진 인코딩.
  • 인간 게놈 WGS BAM 파일 100GB가 CRAM으로 저장하면 30～50GB로 줄어듭니다.
  • 참조가 필요: CRAM 파일을 다시 열려면 원본 참조 fasta가 있어야 합니다.

BAM ↔ CRAM 변환

bash
# BAM → CRAM (참조 지정 필수)
samtools view -T hg38.fa -C -o SAMPLE01.cram SAMPLE01.sorted.bam
samtools index SAMPLE01.cram   # .crai 생성


# CRAM → BAM (다시 되돌리기)
samtools view -T hg38.fa -b -o SAMPLE01.bam SAMPLE01.cram

-T hg38.fa: 참조 지정. 반드시 정렬 때 쓴 그 참조여야 합니다. 그 밖에는 BAM 사용법과 동일합니다.

CRAM 실무 함정 세 가지

  1. 참조 유실 시 재현 불가: CRAM 파일만 남기고 참조를 지우면 리드 원본을 다시는 못 봅니다. 참조도 아카이빙 필수.
  2. 참조 버전 불일치: hg19로 정렬한 뒤 hg38로 열려 하면 각 위치 서열이 안 맞아 오류. 참조 SHA1을 CRAM 헤더에 기록해두는 게 원 규격.
  3. 일부 옛 도구가 CRAM 미지원: 새 파이프라인은 대개 지원. 그래도 참조 소스가 명확한 곳에만 CRAM 저장 원칙.

실무 파이프라인 흐름 — 한 장으로

BWA-MEM에서 GATK 진입 전까지 samtools 사용 순서를 한 장에 정리합니다.

bash
# 1) 정렬 + 좌표 정렬 (스트리밍)
bwa-mem2 mem -t 16 -R "$RG" hg38.fa R1.fq.gz R2.fq.gz \
  | samtools sort -@ 4 -o S1.bam -


# 2) 인덱스
samtools index S1.bam


# 3) 파이프라인 대시보드
samtools flagstat S1.bam > S1.flagstat.txt
samtools idxstats S1.bam > S1.idxstats.tsv


# 4) 아카이빙용 CRAM 변환
samtools view -T hg38.fa -C -o S1.cram S1.bam
samtools index S1.cram


# 5) 파이프라인 진입: MarkDuplicates(S07) 는 sorted BAM으로 넘김

30개 샘플이라도 이 다섯 단계가 뼈대입니다.

손 계산 — CRAM이 왜 그렇게 작아지나

인간 WGS 커버리지 30배 리드가 3억 개, 각 리드 150bp라 하면 BAM에 인코딩된 서열만 45GB입니다. 그런데 리드 하나당 참조와 다른 염기는 평균 몇 개일까요.

  • 인간 개인 게놈의 다형성: 약 4,000,000 SNP (게놈 30억 bp의 0.13%).
  • 리드 150bp에 SNP는 통계적으로 0.2개.
  • 즉 리드 하나에 참조와 다른 위치가 1개 미만입니다.

CRAM은 이 "1개 미만"만 저장하므로 서열 정보량이 원래의 1% 미만으로 줄어듭니다. 나머지 90% 감소분은 헤더 인코딩과 품질 문자열 손실 압축(옵션)에서 옵니다. 저장 공간을 아끼면서도 리드를 다시 복원할 수 있다는 것이 CRAM 설계의 아름다움입니다.

Colab에서 손을 움직여봅시다

bash
!apt-get install -y samtools >/dev/null
!wget -q https://sra-pub-src-1.s3.amazonaws.com/SRR1039508/SRR1039508_1.fastq.gz -O r1.fq.gz
bash
# 미리 정렬된 소형 BAM으로 명령 감각 잡기
!wget -q https://github.com/samtools/samtools/raw/master/examples/toy.bam
!samtools view -h toy.bam | head
!samtools flagstat toy.bam
!samtools index toy.bam
!samtools view toy.bam ref:100-500

이 소형 예제 하나로 view · flagstat · index · 구간 쿼리까지 5분 안에 감을 잡을 수 있습니다.

자주 쓰는 원라이너 다섯 개

파이프라인 로그에서 자주 조회하는 명령들입니다.

bash
# 특정 유전자 영역만 뽑아 새 BAM으로
samtools view -h -b in.bam chr17:41196312-41277500 > BRCA1.bam


# unmapped 리드만 FASTQ로 뽑기 (오염·크로스 종 확인용)
samtools view -b -f 4 in.bam | samtools fastq - > unmapped.fq


# duplicate 표시 리드만 세기
samtools view -c -f 1024 in.bam


# 리드 그룹별 커버리지 (샘플별 요약)
samtools view -H in.bam | grep '^@RG'


# 두 BAM 병합
samtools merge -@ 4 merged.bam S1.bam S2.bam S3.bam

이 다섯 개는 매주 씁니다.

CS 매핑

  • 레퍼런스 기반 압축: 원본 대신 차이만 저장. DryBench의 "레퍼런스 기반 압축" 편(reference-based-compression) 참고.
  • BGZF (Block Gzip Format): BAM은 gzip 블록을 잘라서 저장 → 임의 위치 접근 가능. 인덱스는 이 블록 오프셋을 담습니다.
  • 스트리밍 CLI: samtools의 -(stdin/stdout)로 파이프 연결 → Unix 철학의 정수.

마무리

이 다섯 명령이면 BAM 세계의 95%가 조작 가능합니다. 다음 편(S06)에서는 GATK4 Best Practices의 큰 그림을 잡고, S07에서 MarkDuplicates로 지금까지의 sorted BAM을 이어받습니다. samtools는 GATK 파이프라인에서도 앞뒤로 계속 얼굴을 내밉니다.

더 깊게 파고 싶다면

  • Danecek P. et al. (2021), Twelve years of SAMtools and BCFtools. GigaScience 10:giab008. samtools 12년 개발사와 설계 결정 총정리.
  • Hsi-Yang Fritz M. et al. (2011), Efficient storage of high throughput DNA sequencing data using reference-based compression. Genome Research 21:734. CRAM의 원 논문.
  • htslib 공식 문서 (https://www.htslib.org/): API·CLI 규격 원전.
  • EMBL-EBI Training — Working with BAM files 무료 튜토리얼. 각 명령의 use case 예시가 풍부합니다.
  • Broad Institute BroadE — Working with SAM files 유튜브. samtools view의 실제 리드 조사 사례.

BAM 하나 열어 flagstat과 view를 눌러본 사람만이 다음 GATK 편의 옵션에 자신감을 가집니다. Colab에서 5분만 손을 움직여봅시다.