2022년 2월 21일 월요일

FASTA file을 탭으로 구분된 파일로 펼치기(fasta2tbl 스크립트 + 줄 번호 삽입)

그렇게 중요하지 않은 일에 과도하게 집착을 하는 것은 아닌지 모르겠다. 어제 올린 글('목록을 이용하여 FASTA file에서 원하는 서열을 뽑아내기')에서 이미 awk를 이용한 방법을 소개하였었다. 이것은 특별한 문제가 없이 잘 작동하지만, description 필드를 구성하는 과정에서 사소한 오류가 있다. Description을 포함하는 서열 ID 라인을 처리할 때, awk의 필드 변수인 $1을 null 문자로 바꾸어 버린다. 그러고 나서 $0을 출력하게 만들게 하였더니 $1만 제거되고 그 뒤를 따르는 공백 문자가 description의 시작 부분에 남아 보기가 싫다. 때로는 description 없이 서열 ID만 존재하는 FASTA file도 있을 것이다.

많은 부분을 차지하지 않는 예외 사항까지 잘 처리할 수 있는 스크립트를 짜는데 의외로 많은 시간이 걸린다. 그리고 나는 이러한 마무리 과정에서 오히려 큰 희열을 느낀다. 이것이 바람직한 코딩의 자세인지는 잘 모르겠다. 보람을 상회하는 희열 추구는 자가발전의 경향이 있어서 에너지를 많이 소모하게 만드는데...

인터넷을 뒤져서 가장 완성도가 높아 보이는 awk 스크립트를 구한 다음, 기능을 확장하였다. 아래의 것을 fasta2tbl이라는 파일로 저장한 뒤 처리할 FASTA file 이름을 인수로 주어서 실행하면 된다. 결과물은 서열 ID, description(없는 경우 'NO_DESC'로 표시), 그리고 서열로 이루어진 3개의 컬럼으로 표현된다.

#!/usr/bin/awk -f
# FASTA file을 seq_ID, description, sequence의 세 컬럼으로 이루어진 tab-delimted file로 전환한다.
# description이 없는 경우 두 번째 컬럼에 'NO_DESC'을 채워 넣는다.
{
        if ($1~/^>/) {
                id = $1;
                sub(/^>/, "", id);
                if (NF == 1 )
                    desc = "NO_DESC"
                else
                    desc = $0;
                    sub(/^>[^\s]+\s/, "", desc)
                if (NR>1)
                        printf("\n%s\t%s\t", id, desc)
                else
                        printf("%s\t%s\t", id, desc)
        } else {
                printf("%s", $0)
        }
} END { printf "\n" }

# adapted from Josep Abril's Fasta2Tbl
# https://bioinformatics.stackexchange.com/questions/2649/how-to-convert-fasta-file-to-tab-delimited-file

awk는 세미콜론을 쓰는 방식이나 중괄호의 사용 등이 그렇게 엄격하지는 않은 것 같다. 문법이 헐렁하면 사용자는 게을러지는데...

다음 명령어 한 줄이면 FASTA file을 TSV file로 전환한 뒤 라인 번호(with leading zeros)를 삽입한 네 컬럼짜리 결과 파일을 얻게 될 것이다.

$ fasta2tbl input.fa | nl -n rz - > SeqWithLineNum.tab

댓글 없음: