Fork me on GitHub 단맛만좋아요 :: 단맛만좋아요 (2 Page)

그래프를 이용한 태스크 표현

흔히 바이오인포매틱스 분석이라고 하는 경우 스크립트나 모듈을 작성하여 일련의 분석을 수행하곤 한다. 그러나 단순한 형태의 일이 아니라 더욱 고난이도의 일을 처리하다가 보면 (물론 대부분이 그렇지만) 태스크의 의존성을 고려해야 하는 경우가 많다. 그래프를 이용하여 이를 표현해 보면서 어떻게 의존성과 결함을 고려한 스케줄러를 만들 수 있는지 생각해 보도록 하자.



총 4개의 노드 (Task1, Task2, Task3, Task4)와 엣지로 구성된 그래프로 각각의 엣지는 다음과 같은 의존성을 가진다.


Task1=>[Task2, Task3]

Task2=>[Task4]

Task3=>[Task4]

Task4=>[]


즉, Task1이 끝나야 Task2,3이 수행되고 Task4는 Task2,3가 끝나야 수행된다. 이렇게 각각의 노드(태스크)는 의존성(엣지)을 가지고 최종적으로 위의 그래프는 Task1-Task2, Task1-Task3, Task2-Task4, Task3-Task4로 표현될 수 있다.

워크플로우 상호 의존 태스크 관리

SGE를 비롯한 대부분의 스케줄러는 workflow interdependence (워크플로우 상호의존) 개념을 가지고 있으며 다음과 같이 표현될 수 있다.

qsub -N Task1
qsub -N Task2 -hold_jid Task1
qsub -N Task3 -hold_jid Task1
qsub -N Task4 -hold_jid Task2, Task3


어때유 간단하쥬,

결함을 허용하는 워크플로우 만들기

그런데 문제는 어떻게 결함을 허용하게 하느냐는 것이다. 여러가지 방법이 있겠지만, 필자는 각 Task의 결과 파일과 결과 파일의 MD5 체크섬을 이용하는 방안에 대해서 이야기 해 볼 것이다.


각 Task는 Task의 결과로 생성되(어야하)는 ouput file이 존재한다. 따라서 Task의 실행전에 해당 output file이 존재하는 또한 저장된 ouput file의 MD5 값이 outputfile을 실제 MD5를 실행해서 얻은 값과 동일한지(결과값이 변경되었는지를 확인하기 위해) 를 확인한다. 만약 두 조건을 만족한다면 해당 Task는 이미 정상적으로 실행되었던 것이기에  skip하고 그렇지 않은 경우에는 해당 Task가 처음 수행되거나 또는 비정상적으로 종료된 것이기 때문에 Task 고유의 action을 수행한다.


if [[ -f $OUTFILE && -f $MD5FILE && `md5sum ${OUTFILE} | cut -f1 -d" "` = `cat ${MD5FILE} | cut -f1 -d" "` ]]; then

   skip

else

   do task action

   task_md5=`md5sum $OUTPUTFILE > $MD5FILE

fi


위 그림과 같이 Task3이 어떤 이유에서든 실패했다면, 실패의 원인을 찾아 해결한 후 다시 모든 Task job을 resubmit하면 Task1은 skip이 되고 (이미 정상적으로 수행되었기때문) Task3부터는 재실행하게 되어 결국 자원낭비 없이 모든 Task들이 수행된다.

Job Array를 이용한 반복 태스크 관리

한가지 더 이쪽 분야의 일을 하다가 보면 특정한 일을 여러번 반복해야하는 경우가 생긴다. 가령 염색체별로 23번을 반복해야한다거나 샘플별로 동일한 작업을 여러번 해야 한다거나 말이다. Task2, Task3...TaskN의 작업은 변수1개만 다를뿐 동일한 경우에는 Job Array를 이용한다.


qsub -N Task2toN -t 1-100 


-t 옵션을 주어 1부터 100까지 100번의 job을 차례로 submit한다. 그럼 나는  job 스크립트에 다음과 같이 array별로 변수를 바꾸어 할당하여 job을 수행하도록 하면 된다.


sample=awk "NR==${SGE_TASK_ID}" samplelist.txt

java -jar GATK.jar -T GenotypeGvcf -I ${sample}.bam -O ${sample}.vcf


samplelist.txt에는 sample 이름 100개가 담겨져 있고 스크립트는 samplelist.txt 파일을 한줄한줄 읽어 두번째  모든 샘플에 대해 GATK 작업을 100번 수행한다. 100번을 job을 submit하지 않고도 100번을 job을 수행할 수 있는 것이다. 이렇게 간단하게 embarrassingly parallel를 손쉽게 구현이 가능하다.


간단하게 바이오인포매틱스 워크플로우에서 SGE 스케줄러를 이용하여 의존성 있는 Task들에 대해서 결함을 허용하는 한편 간단히 array job을 수행하는 방법에 대해서 알아보았다.


실제로는 NoFlo 등을 이용하여 프론트엔드단을 구현하고 그안에 SGE 스케줄러와 위에 설명한 방법들을 결합하여 워크플로우를 작성하여 실행하도록 하면 끝.

Job 모니터링과 fail job 복구

마지막으로 스케줄러를 사용하다보면 간혹 shared storage가 네트워크의 과부하 등으로 제대로 작동되지 않는 경우 SGE는 job의 상태를 Eqw로 빠드리고 에러를 확인해 보면 "can't chdir to directory: No such file or directory"라고 shared storage에 접근하지 못한다는 에러를 확인 할 수 있다. 이때는 아까도 말했듯이 간헐적인 네트워크 과부하 등의 원인으로 shared staroge 접근에 임시적으로 문제가 생길 수 있으니 이때는 job을 다시 실행하도록 하는데 qmod의 -cj 옵션으로 에러를 해결(clear)했으니 다시 시도하라는 것이다.


간단히 다음의 bash 스크립트로 job 상태를 xml 형태로 받아 파싱하여 queue의 상태를 모니터링하도록 하자.


로그 관리

분석을 수행하면서 생기는 로그는 다음의 3가지로 나누어 볼 수 있다. 스케줄러에서 생성되는 Job별 STDOUT, STDERR과 전체 워크플로우에 사용자가 넣은 로그로 각 로그는 웹 상에서 자동으로 갱신되어 보여지도록 한다. 이때는 js-logtail 등을 이용하면 쉽게 로그 파일을 처리할 수 있다.





저작자 표시 비영리 동일 조건 변경 허락
신고
Software enginner of GenomeCloud. Covers bioinformatics, computational biology, and life science informatics.
Posted in : 빅데이터분석 at 2015.10.13 18:58
Currently 2 comments want to say something now?
BGZF (Blocked GNU Zip Format)
2015.07.01 16:18 | 바이오인포매틱스

Random Access

BAM 파일의 경우에는 BGZF를 이용하기 때문에 원하는 곳으로 빠르게 access가 가능하다. 우리가 흔히 사용하는 GZIP (GNU ZIP) 보다는 압축효율 (압축했을때 용량)이 떨어지지만 random access가 가능하다는 잇점으로 인해 BAM 파일(BAM의 경우 재빠르게 자신이 원하는 position을 뷰잉하는데 많이 사용하기 때문)에서 사용하는 기술이다.


용량이 큰 텍스트 파일을 압축해 놓고 파일의 어느 부분이던지 랜덤하게 액세스 가능하기 때문에 그 활용도가 높은데 특히나 클러스터를 이용하는 경우 파일을 분할하는 등의 I/O 작업이 필요 없기 때문에 그 활용이 매우 높다고 할 수 있다. 여러 활용중 하나로 FASTQ 파일에 적용하여 사용하고 있다. 

FASTQ 응용

일루미나는 Base Call을 수행한 후 기본적으로 GZIP을 이용하여 압축하게 되는데 이를 다음과 같은 명령으로 BGZF로 변환한다. 뭐 언뜻 보기에는 BGZF가 GZIP과 호환되기 때문에 이게 BGZF인지 확인해 봐야 한다. grabix 툴의 check를 이용하면 간단하게 알 수 있다.


# zcat TN1506L0327_2.fq.gz |bgzip /dev/stdin -c > TN1506L0327_2.fastq.gz

# /BIO/app/grabix/grabix check TN1506L0327_2.fastq.gz

yes

# /BIO/app/grabix/grabix index TN1506L0327_2.fastq.gz

# /BIO/app/grabix/grabix grab TN1506L0327_2.fastq.gz 1 10

@HWI-ST840:481:HVKTWADXX:1:1101:1213:1943 2:N:0:ATTGGCTC

GGTACCAGTACCATGCTGTTTTGGTTACTGTAGCCTTGTAGCACAGNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN

+

BB@DFFFFFHHHHJJJJJIJJJJJIIJJJJIIJIJJJJJJJJIJJI#######################################################

@HWI-ST840:481:HVKTWADXX:1:1101:1986:1934 2:N:0:ATTGGCTC

ATATCTACACCACGCAAAGTGATGTGTAAGTGTGGGTGTTGCTCNCNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN

+

@@CDDFFAHFFHHIIIIIIHHIHECDHIIEFHHHGI@DDGIIIG#########################################################

@HWI-ST840:481:HVKTWADXX:1:1101:1874:1947 2:N:0:ATTGGCTC

CCTAGGAATCCAACTTACAAGGGATGCGAAGGACCTCTTCAAGAAGNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN

BWA 활용

이렇게 BGZF 포맷의 FASTQ 파일은 굳이 압축을 풀지 않아도 원하는 라인을 바로 확인이 가능하게 된다. 이렇게 되면 우리는 클러스터 노드에 다음과 같이 직접 파일을 분할 할 필요 없이 "grabix grab 파일명 start_line end_line"을 통해 다수의 노드에 BWA 매핑 작업을 분할하여 수행할 수 있다.


# bwa mem -M -t 7  -R '@RG\tID:1\tPL:illumina\tPU:1_2014-08-13_dream-syn3\tSM:syn3-normal' -v 1 GRCh37.fa <(grabix grab normal-1.fq.gz 20000001 40000000) <(grabix grab normal-2.fq.gz 20000001 40000000)


많은 생물학 데이터 FASTQ, FASTA, SAM, BED, XML, Tab 구분 파일 등등 모두 적용이 가능하며 다양하게 활용이 가능하다. 물론 여러개의 BGZF 파일을 압축해제 과정이 필요없이 하나의 BGZF로 묶어주는 방법도 장혜식 박사의 블로그에 있으니 자신의 데이터에 적용해서 활용해 보기 바란다.

Python을 이용한 BGZF 헤더 다루기

RFC1952 문서에 정의된 GZIP에서 extra field를 활용하는 BGZF는 extra field의 두개의 subfield id에 각각 66, 67의 값을 가지며, EOF 마커로 "1f 8b 08 04 00 00 00 00 00 ff 06 00 42 43 02 00 1b 00 03 00 00 00 00 00 00 00 00 00"의 16진수를 가진다.

 


Python을 이용하여 파일을 읽고 BGZF 포맷이 맞는지 확인하는 코드를 작성해 보도록 하자. 대학원 시절 프로토콜 정의하고 구현하던게 갑자기 생각이 난다. 암튼 뭐 어쨌건간에 우선 파일의 첫번째 블록의 헤더 부분을 읽어오는데 기본 헤더와 함께 SI1,SI2의 extra subfield의 일부분을 읽는다. 읽어오는 부분은 총 14바이트 하나씩 출력해보면 위에 정의된 값들이 출력되는 것을 확인할 수 있으며, XLEN 부분은 총 6바이트가 extra subfield로 지정된 것을 확인 할 수 있다.


from struct import unpack


bgzf_file=open('TN1506L0327_1.fastq.gz')

header = bgzf_file.read(14)


id1, id2, cm, flg, time, xfl, os,extra_len, extra_id1, extra_id2 = unpack('<BBBBLBBHBB', header)


>>> from struct import unpack

>>>

>>> bgzf_file=open('TN1506L0327_1.fastq.gz')

>>> header = bgzf_file.read(14)

>>>

>>> id1, id2, cm, flg, time, xfl, os,extra_len, extra_id1, extra_id2 = unpack('<BBBBLBBHBB', header)

>>> id1

31

>>> id2

139

>>> cm

8

>>> flg

4

>>> extra_len

6

>>> extra_id1

66

>>> extra_id2

67


BGZF 파일의 EOF를 읽어들여 파일이 제대로 생성된 것인지를 확인해 본다. 이전에 언급했듯이 파일의 마지막 28바이트를 읽어들여 EOF 마커와 동일한지 확인한다. 왜 BGZF의 첫블록을 읽고 마지막 EOF를 확인하는지는 많은 데이터가 짤리는 경우가 많아서 확인 과정을 거친 후 진행하는게 건강에 좋아서랄까? ㅎ


>>> bgzf_file.seek(-28,2)

>>> eof=bgzf_file.read(28)

>> eof

'\x1f\x8b\x08\x04\x00\x00\x00\x00\x00\xff\x06\x00BC\x02\x00\x1b\x00\x03\x00\x00\x00\x00\x00\x00\x00\x00\x00'


저작자 표시 비영리 동일 조건 변경 허락
신고
Software enginner of GenomeCloud. Covers bioinformatics, computational biology, and life science informatics.
Posted in : 바이오인포매틱스 at 2015.07.01 16:18
Currently 댓글이 없습니다. comments want to say something now?

티스토리 툴바