NGS 데이터 분석 - GATK 전체 파이프라인 및 설치

GATK (Genome Analysis Toolkit)

Broad Institute에서 개발한 것으로, 인간 유전체와 엑솜 데이터를 중심으로 한 variant discovery와 관련된 다양한 기능을 제공한다. Germline DNA의 SNPs과 Indels를 발굴을 시작으로 somatic short variant calling, CNV, SV까지 확장하며, Processing이나 quality control와 같은 다양한 기능을 제공한다.

핵심은 Sequencing data에서 Genomic variation을 찾고자 하는 것

GATK Best Practices

GATK에서는 End-to-end Reads-to-Variants 워크플로우(link)를 제공하고 있으며, 버전 4.0부터 GATK 내에서 Picard를 모두 사용할 수 있다. 개인 노트북이나 서버에 설치하거나 Broad Institute에서 제공하는 클라우드 플랫폼인 Terra를 통해 Jupyter Notebook 상에서 분석을 할 수 있다.

GATK Best Practices를 통해 raw 시퀀싱으로부터 분석에 적합한 filtered variant callset을 생성하며, 실험 설계나 QC, 파이프라인 옵션 등에 대해 가이드라인을 제공한다. 따라서 각자의 실험 설계에 알맞는 적절한 파이프라인을 선택해야 한다 (예: 플랫폼, DNA- or RNA-Seq, WES or WGS (PCR-free or PCR+), paired- or single-end, read length, expected average coverage, somatic data 등). 이때 기본적으로 Illumina를 기반으로 개발되었다는 점과 다양한 요소(시퀀싱 기술, 하드웨어 인프라 등)가 영향을 미친다는 점을 반드시 고려해야 한다.

0. 분석 단계

상황에 따라 다음과 같이 2-3개의 단계로 나눌 수 있다.

(1) 데이터 전처리 (Data Pre-processing) : 모든 분석의 첫 단계로, raw 시퀀싱 파일(FASTQ, uBAM)로부터 BAM 파일을 생성한다. Reference 유전체에 대해 alignment를 하고, 기술적인 biases를 고려하여 분석에 적합한 형태로 데이터 클리닝하는 과정이 진행된다.

(2) 변이 발굴 (Variant Discovery) : BAM 파일로부터 variant calling을 실시한다. 개인으로부터 유전적 변이를 발굴하고 실험 설계에 맞은 필터링 방법을 적용한다. 기본적으로 VCF 파일이 생성되지만, CNV와 같은 변이의 경우 다른 형태로 결과물을 저장한다.

(3) 추가 분석 : 필터링이나 annotation이 추가적으로 진행될 수 있으며, 이미 알려진 변이 정보를 활용하거나 메타데이터를 통해서 정확도를 평가하고 개선한다.

1. 데이터 전처리

개별 샘플에 대하여 각각 전처리를 시행한다. uBAM 파일을 input으로 하며, FASTQ 파일의 경우 uBAM으로 변환이 가능하다.

데이터 전처리는 다음과 같이 3가지 과정으로 나눌 수 있다.

(1) Map to Reference : 먼저 reads를 reference 유전체에 mapping한 SAM/BAM 파일을 생성하고, coordinate에 따라 정렬한다. (BWA, MergeBamAlignments)

(2) Mark Duplicates : PCR amplication bias 등을 완화하기 위하여 duplicates를 찾는다. (MarkDuplicatesSpark / MarkDuplicates + SortSam)

(3) Base Quality Scores Recalibration : Variant calling 알고리즘은 개별 base calls에 부여된 quality score에 크게 의존하므로, Base quality scores를 recalibrate한다. (BaseRecalibrator, ApplyBQSR, AnalyzeCovariates (optional))

 

[참고]

• GATK 시작하기 (link)
• GATK Best Practices에 대하여 (link)
• 데이터 전처리 (link)

2. 변이 발굴

연구 코호트를 구성하는 각 샘플별 BAM 파일을 GATK을 통해 데이터 전처리를 시행한 후 다음 변이 발굴을 진행한다.

Variant Discovery는 Short variant (SNPs + Indels)와 Copy number variant (CNVs) 두가지로 나눌 수 있다.

시퀀싱의 경우 DNA-seq (Germline, Somatic, Mitochondrial)와 RNA-seq이 있으며, 세부적으로 WGS or WES, gene-panels 등 다양한 실험 설계가 존재한다.

A. Germline (short variant, copy number variant)

B. Somatic (short variant, copy number variant)

C. RNA-seq (short variant)

D. Mitochondrial (short variant)

Germline Short Variant Discovery

A. Cohort Data

Germline 코호트 데이터인 경우,

(1) Call Variants Per-Sample : 샘플별로 variant calling하여 GVCF 포맷을 생성한다 (HaplotypeCaller)

(2) Consolidate GVCFs : 그 후 여러 샘플의 GVCFs을 GenomicsDB로 합친다 (GenomicsDBImport)

(3) Joint-Call Cohort : Joint genotyping을 시행한다 (GenotypeGVCFs)

(4) Variant Quality Score Recalibration & Filter Variants : 적절한 정밀도와 민감도에 따라 VQSR 필터링을 통해 마지막 multisample callset을 얻는다 (VariantRecalibrator, ApplyVQSR)

 

이때, 동시에 많은 샘플이 주어질수록 Germline variant discovery 알고리즘의 정확도와 민감도가 매우 향상된다. 특히 코호트의 모든 샘플에 걸쳐서 잠재적인 변이 위치를 나타내는 유전형 행렬을 구축하는 것이 필요하며, 이는 샘플별로 따로 생성된 variant calling를 합치는 기존의 원시적인 접근과는 다르다.

초기 variant discovery 단계에서는 multiple per-sample BAM 파일을 joint anlaysis 시행하였지만, 샘플수에 따라 조정이 잘 되지 않아 연구 코호트의 크기에 매우 제한을 받게 되었다. 게다가 새로운 샘플이 추가될 때마다 variant calling을 다시 실시해야한다는 점에서 점진적으로 샘플수를 늘려나갈 수 없었다. 따라서 버전 3부터 새로운 접근 방식이 도입되어 기존의 variant calling을 2개의 내부 과정으로 나누게 되었다. 이를 통해 점진적으로 샘플수를 늘려 대규모 코호트에서도 분석이 가능하도록 하였다.

(1) 샘플별 variant context 통계와 모든 genotype likelihoods의 계산 : BAM 파일이 필요하며 computationally expensive

(2) 샘플별 genotype posterior probabilities 계산 : 코호트의 모든 샘플의 genotype likelihoods를 기반으로 하며 computationally cheap

B. Single Sample Data

Germline 단일샘플 데이터인 경우,

(1) Call Variants Per-Sample : 단일 BAM 파일을 variants calling을 하여 VCF 파일을 생성한다 (HaplotypeCaller; 기본 단일샘플 모드)

(2) Filter Variants : 모델을 통해 각 변이의 quality 예측값으로 annotate하고, 해당 값을 기반으로 filtering을 실시한다 (CNNScoreVariants, FIlterVariantTranches)

Somatic Short Variant Discovery

Somatic인 경우, 대상자의 tumor & normal 샘플의 BAMs(T-N BAMs)을 통해 tumor 샘플에서의 variant를 찾고자 한다. 이를 위해 먼저 대규모로 somatic variants 후보를 만든 후, filtering을 실시하여 confident set을 얻는다.

(1) Call candidate variants (Mutect2)

(2) Calculate Contamination (GetPileupSummaries, CalculateContamination)

(3) Learn Orientation Bias Artifacts (LearnReadOrientationModel)

(4) Filter Variants (FilterMutectCalls)

(5) Annotate Variants (Funcotator)

Somatic Copy Number Variant Discovery

환자 샘플에서 somatic CNVs 발굴을 위해서는 적절한 Panel of Normals (PON)이 필요하다. PON은 Normal BAMs를 통해 만들 수 있다.

RNA-seq Short Variant Discovery

RNAseq일 경우, 여러 uBAM 파일을 한번에 분석하며 joint calling은 지원하지 않는다.

(1) Mapping to the Reference (STAR)

(2) Data Cleanup (MergeBamAlignment, MarkDuplicates)

(3) SplitNCigarReads (SplitNCigarReads)

(4) Base Quality Recalibration (BaseRecalibrator, ApplyBQSR, AnalyzeCovariates)

(5) Variant Calling (HaplotypeCaller)

(6) Variant Filtering (VariantFiltration)

 

---

GATK 설치하기 (Docker)

• (How to) Run GATK in a Docker container (link)
• Dockerhub (link)
• Github : Official code repository for GATK versions 4 and up (link)

리눅스 기반으로 실행되며 MacOS에서 사용이 가능하지만 윈도우는 지원하지 않는다. Docker을 사용하기를 권장하고 있으며, Java 1.8가 필수이며, R이나 Python이 사용되기도 한다.

Docker 설치 후, Linux shell에서 Docker hub로부터 gatk 이미지를 다운받은 후, 공유폴더를 지정하여 container을 생성한다.

docker --version # Docker version 20.10.7
docker pull broadinstitute/gatk:4.2.4.1
docker run -v ~/my_project:/gatk/my_data -it broadinstitute/gatk

Container 내부에서 gatk 명령어를 통해 원하는 tool과 argument를 지정하여 프로그램을 돌린다.

./gatk --help # For help on using gatk itself
./gatk --list # Print a list of available tools
./gatk ToolName --help # Print help for a particular tool
./gatk ToolName toolArguments # Run a non-Spark tool, or to run a Spark tool locally

# Example
./gatk PrintReads -I input.bam -O output.bam
./gatk PrintReadsSpark -I input.bam -O output.bam