유전자정보분석

Genome Analysis with MapReduce

hongiiv 2012. 8. 23. 21:46
반응형
여: 옵빠 나 밤(BAM)가지고 있는데 요것들을 합치려 하니깐 리드들의 아이디가 서로 겹치는게 있어서 에러나 이거 고치는데 시간 오래 걸려?
남: 음! 그것두 졸라 한 오백 라인정도 하루죙일 코딩해야해! 
(샘포맷으로 바꾸고 헤더 정보 읽은 담에 리드 그룹을 파싱해서 아이디로 해쉬 테이블 만들어 넣고 리드들을 루푸 돌아가면서 해쉬에 있는 리드 그룹 아이디를 찾고 플랫폼 태그를 찾아서 리드 이름에 추가하고 다시 밤포맷으로 바꾸고 인덱싱 새로하고...)

"The presentation from Eli Lilly is a great introduction to developing your own custom GATK Walkers in Java."

GATK Overview
자! 당신이 NGS 데이터를 다루는데 좋은 아이디어를 GATK를 이용해서 구현하는지에 대한 것이 바로 오늘 준비한 이야기이다. 그럼 기존에 어떠한 상황에서 GATK를 사용했었는지에 대해서 우선 짚고 넘어가보도록 하겠다. 다시 한번 Genome Analysis Toolkit (GATK)는 Broad Institute에서 개발된 Application Programming Interfae (API)와 NGS 데이터를 분석하기 위한 소프트웨어 툴의 집합이다.GATK는 alignment, FASTA 파일 통계 분석, 샘플 커버리지 계산등의 서로 다른 컴포넌트로 구성되어 있다. 

GATK의 Variant 발굴은 4개의 단계로 나누어 볼 수 있다.

1. Initial Mapping
2. Read Refinement/Recalibration
3. SNP Calling
4. SNP Filtering/Variant Quality Score Recalibration

Read의 Refinement와 Recalibration을 위해서 GATK의 CountCovariatesTableRecalibration 툴이 사용되며, SNP Calling에는 UnifiedGenotyper, SNP calling 후에는 VariantRecalibratorApplyRecalibration이 사용되며, 위의 단계에는 없지만 annotation을 위해서 VariantAnnotator를 각각 사용한다. 아마도 GATK는 Human에 있어서 현존하는 최고/최신의 Variant Calling 툴이 아닐까한다.

임상으로서의 NGS 데이터
그럼 간단히 임상에 NGS 데이터가 활용되기 위해서 필요한 조건들이 무엇인지 한번 생각해보자? highly accurate clinical grade sequencing and variation detection method!

Clinical grade genome
  • 98 percent genome coverage
  • 1 error per million base (SNPs + small indels)
  • Full haplotype phasing
  • Structural variations
Genome Analysis with MapReduce
위에서 언급한 임상에 적용하기 위해 데이터를 분석하는데에 있어서 대용량의 NGS 데이터를 다루는데에는 GATK가 유용하다는 것을 알았다. 그러면 이제 좀 더 자세히 알아보도록 하겠다.

게놈시퀀싱 기술의 발달로 NGS를 통해 생성되는 데이터는 기하 급수적인 속도로 늘어났다. 일례로 NGS를 이용한 1000 Genomes 프로젝트가 있다. 여기서는 앞서 누누히 이야기했듯이 GATK라는 구조화된 프로그래밍 프레임워크를 이용하여 NGS 데이터를 분석하는 것에 대해서 알아볼 것이다.

GATK and MapReduce
이 프레임워크는 MapReduce라는 프로그래밍 패러다임을 사용하여 쉽고 효율적으로 게놈 데이터를 분석하기 위한 도구를 만드는데에 사용된다. MapReduce는 대용량의 데이터를 처리하는데에 사용되는 분산 프로그래밍 프레임워크이다. 

클러스터 (or 클라우드) 컴퓨팅상에서 MapReudce는 대용량 데이터를 저장하기 위해 Hadoop 분산 파일 시스템을 제공한다. Hadoop을 이용한 간단한 WordCount와 Hadoop-BAM이라는 응용 프로그램을 통해 Hadoop과 MapReduce를 사용하는 방법과 GATK에 있는 간단한 예제를 통해 어떻게 MapReduce를 구현하는지에 대해서 설명할 것이다.

병렬화는 GATK의 특징GATK는 human medical sequencing project (Cancer Genome Atlas)에 사용되었는데, NGS 데이터를 분석하고 분석 알고리즘에서 데이터 액세스에 관련한 패턴을 분리하여 데이터 관리 문제를 해결하는데 도움을 준다. (뒤에 설명하겠지만 traversals와 walkers)

GATK의 특징 중 하나는 레퍼런스 데이터와 snp에 대한 정보를 포함하는 엄청난 시퀀스 데이터를 분할 한다는 것이다. 또 다른 특징은 사용자가 병렬처리를 할 수 있도록 옵션을 제공한다는 것이다. 즉 분할하여 처리가 가능하기 때문에 병렬처리가 가능해진다. 이뿐만 아니라 몇개의 BAM (즉 여러 샘플) 파일이나 시퀀스 파일을 한번에 분석 할 수 있다는 것이다. 

GATK의 병렬화 방법
이러한 병렬화는 공유메모리상에서 자동으로 병렬화하는 방법과 하나의 machine에서 다수의 traversal 엔진과 walker과 다수의 인스턴스를 통해 가능(Scatter/Gather)하다.

GATK 아키텍처
GATK는 traversals와 walkers로 구성되는데 traversals는 walkers가 접근하기 위한 데이터에 대한 정보를 제공하고 walker는 map과 reduce를 함수를 이용하여 데이터를 다루는 역할을 한다.

read-based traversals/ReadWalker
시퀀스 read와 이와 연관된 reference 데이터에 대해서 반복하면서 접근한다. walker는 한번에 하나의 read를 처리하게 된다. 이러한 read baesd traversals과 ReadWalker를 이용하여 다음과 같은 작업들을 수행할 수 있다. (1) read에 들어있는 metadata에 접근하여 이를 셋팅 (2) indel을 찾거나 realign을 수행

관련된 GATK 프로그램으로는 CycelQualityWalker, TableRecalibrationWalker, CountReadsWalker, IndelRealigner 등이 있으며, 앞서 남녀간의 대화에 나온 내용도 바로 read-base로 처리가 가능하다.



locus-based traversals/LocusWalker
locus 기반은 레퍼런스의 순서에 따른 데이터 접근으로 주어진 locus의 base들에 접근할 수 있다. pipeup 등에 사용될 수 있으며, 하나의 위치상에 있는 레퍼런스, read 들에 접근하게 된다. Variant calling이나 Depth of coverage 계산, 주어진 영역내의 GC content, read error rates 등을 계산하는데 사용될 수 있다. 

Rod-based traversals/RodWalker
locus-based와 비슷하게 Rod-base가 있는데, locus base가 주어진 locus로부터 ref에 따라 순차적으로 접근하것에 비해 vcf 파일과 같이 순차적이지 않고 산발적인 데이터에 접근할 수 있다. 사용되는 것은 locus-based와 유사하며 GATK의 VariantEval, PhaseByTransmission, VariantAnnotator, VariantRecalibrator, SelectVariants등이 있다.



전체적인 GATK의 구조는 GATK가 가장 상단에 위치하고 Traversal engine에 의해  genome 데이터의 특징에 따라 접근을 하고 접근된 데이터를 활용하는 Analysis tool이 존재하게 된다. analysis tool은 이미 GATK에 의해 만들어진 것과 사용자가 직접 구현할 수도 있다.



GATK를 이용한 사용자 Analysis tool 만들기 (1) - ReadWalker
(1) GATK는 소스코드와 미리 빌드(소스코드를 컴파일해서 바로 사용 가능한 버전)해 놓은 두가지 버전을 제공한다. 우리는 소스코드를 수정할 것이기 때문에 소스코드를 git에서 다운로드 한다.
(2) 다운로드한 소스코드는 Intellij IDEA를 이용하여 프로젝트를 생성한다. (사용자에 따라서 Eclipse를 써도 무방)

아래 그림은 GATK의 소스코드를 Intellij로 불러온 모습으로 org.broadinstitute.sting 패키지가 보이고 크게 1) alignment에 관련된 패키지 2) analyzecovariates 분석시 사용하는 covariate관련 패키지 3) GATK의 커맨드라인에 관련된 패키지와 4) 주요 gatk관련 패키지로 나뉜다. gatk 패키지에는 앞서 설명했던 traversals와 walkers 패키지가 있고 좀 재미있게 생겨먹은 phonehome이라는 패키지가 보인다. phonehome은 gatk를 수행하고 난 로그를 저장 관리하는 패키지인데 gatk를 돌리고 나면 관련 로그를 gatk에서 지정한 곳(아마존의 S3 스토리지)에 전송한다. 



gatk 패키지의 walkers에 examples라는 패키지(org.broadinstitute.sting.gatk.walkers.example)를 생성하고 HelloRead라는 클래스를 만든다.

HelloRead 클래스는 ReadWalker를 extend하는 패키지로 map, reduceInit, reduce 메소드를 각각 생성한다. PrintStream을 정의하고 map메소드의 read인자에는 읽혀진 read 정보로 read.getReadName()을 사용하여 읽혀진 read에 대한 이름과 "Hello"라는 문자을 추가햐여 출력한다. 


HelloRead를 추가한 gatk 소스코드를 아래와 같이 "ant dist" 명령어를 통해서 컴파일을 수행한다.



아래 그림처럼 해당 영역 (-L 옵션)에 위치하는 read들에 대해서 Hello라는 문자열을 추가하여 해당 read의 이름을 출력한다.
java -jar /root/broadgsa-gatk-fe71742/dist/GenomeAnalysisTK.jar -T HelloRead -R /RAID/spiral/software/data/reference/hg18/hg18.fa -I /RAID/spiral/BAM/TGP2011D0013_AA.realign.recal.redup.matefix.sorted.bam -L chrM:2204-2205
 
GATK의 실행 명령은 -T 다음에 원하는 GATK의 분석 툴 이름, -R 다음에는 reference, -I 다음에는 분석하기를 원하는 파일 즉 입력파일의 경로, -L에는 특정 영역을 지정할 수 있다 물론 생략하면 입력파일 전체에 대해서 수행한다.

명령을 수행하고 나면 아래와 같은 broad에서 지정한 특정 위치의 아마존 스토리지로 로그 정보가 전달된다. 전달된 정보는 XML 형태로 아래와 보는 바와 같은 다양한 정보를 포함해서 전달된다. GATK가 버전업이 되면서 전달되는 정보가 좀 줄긴 했다.



어떤 버전의 GATK를 사용했는지? 입력데이터는 어떤경로에 어떤 파일인지 등등 다양한 정보가 전달된다. 만약 정보를 전달하지 않으려면 broad에 연락하던가 소스코드를 저처럼 살짝 고쳐서 사용하면 된다. 물론 이는 연구용으로 사용할때만 가능한것이고 상업적인 용도로 사용하려면 broad에 연락은 필수!



여기까지 간단하게 나마 GATK의 구조와 간단한 예제를 하나 살펴보았다. 앞으로 BAM(SAM) 파일이나 VCF 파일을 다룰 일이 많은 분들은 한번쯤 GATK를 사용해서 자신의 원하는 custom tool을 만들어 사용하는 것도 좋을 듯! 다음번에는 좀 더 다양한 예제와 GATK의 병렬화 방법에 대해서 좀 더 자세히 알아보도록 하겠다. 글이 중구난방으로 여기왔다 저기왔다 하는감이 많이 든다. 글을 다듬을 시간도 없고 해선 우선 올리는 것에 만족~~!
반응형