[Genomics data processing] awk를 이용해서 bed format 을 chromosome 별로 분리하기

#basic
awk '{close(f);f=$1}{print > f".bed"}'  [input bed]

#making path is also available
awk '{close(f);f=$1}{print > "testdir/"f".bed"}' ./input_whole_genome.bed

---

이 카테고리에 처음 기록을 하게되어 쓰는 사족:

이런 저런 작업들을 하다보면, 인턴때나 석사 초기에는 직접 간단한 스크립트를 짜서... 데이터를 뽑고 그래도 그 와중에 조금씩 틀리기도 하고 고치기도 하고.. 하면서 작업을 했었던 것 같다. 그러다 조금 짬이 차니 그런 간단한, 반복적인, 언젠가 또 할 법한 작업들의 general한 형태 스크립트를 만들어 홈 디렉토리에 넣어두고 필요할 때 꺼내 조금씩만 수정해서 사용하다가... 또 여러 가지에 익숙해지고 보고듣고 하다보니 아 이거 리눅스 커맨드로도 가능할 것 같은데..? 하다가 찾아지는 꿀팁들을, 사용하고 있던 꿀팁들을 조금씩 정리하려고 한다. (문장이 두서가 없는 것은 기분탓..?) 아무튼, 이렇게 정리해두면 나중에 필요할 때 내 블로그 뒤져서 쓰거나 이렇게 기록하면 또 금방 기억 나겠지! 후후훗 (?)

---

BED format (.bed)은 genome 내에 genomic element의 위치를 표현하는 텍스트 파일 형식으로 탭으로 column을 구분한다. 가장 기본 형태는 chromosome(reference sequence name) - start position - end position 3개 컬럽으로 구성된 형태이며 뒤의 컬럼을 추가로 붙여 방향, 추가 정보 등을 표현할 수 있다. 

  엄밀하게 BED format 상에서 start position은 0-based notation을 사용해 표현하고,(실제 시작 위치 포지션에서 -1 한 위치를 시작 위치로 표현하는 방식, 길이를 계산할 때 유용... 한 것 외에 왜 쓰는지는 잘 모르겠음) 3-column 보다 확장된 형태의 bed format의 경우 몇 번째 컬럼에 방향을 표시하는 등 추가적인 규칙이 있긴 하나 아래의 예시는 그것이 엄격하게 적용되어있지는 않다. 기회가 되면 format 들을 따로 다루는 글을 써 봐야 겠다.

 

#bed format
14  417 417 54.5454545454545    6   5
14  428 428 0   0   10
14  431 431 0   0   12
14  435 435 0   0   12
14  437 437 0   0   13
14  442 442 0   0   14
14  445 445 0   0   14
14  451 451 5.55555555555556    1   17
14  454 454 18.1818181818182    4   18
14  455 455 0   0   18
14  458 458 9.52380952380952    2   19

 

  whole genome 에 대해 만들어진 데이터들은, 특히나 위와 같이 base-level로 만들어진 경우 용량이 엄청나게 큰 경우가 많다. vi로만 열어도 수기가의 메모리를 잡아먹는.. (그럴 땐 사실 vi 를 쓰는걸 지양하라고 배우긴 하지만... ) 이런 데이터를 인풋으로 어떤 프로그램을 실행할 경우, 문제 종류에 따라 전체 프로그램 실행 시간을 단축시키기 위해 병렬화 하기도 하는데, 이런 genomic region 데이터를 나눌 수 있는 가장 좋은 방식은 chromosome 이다. genomic region 데이터 관련 작업들이 대개 한 파일만 필요한 경우 보다는 지역에 대한 annotation이 있는 다른 데이터를 함께 사용하게 된다. 이 때 reference sequence(chromosome)구분 없이 line 수로 나눠버리면 한 파일에 들어있는 reference sequence 모두를 처리하기 위해 전체 리소스를 더 잡아먹는 불상사가 생길 수도 있으므로 보통 chromosome level로 쪼개서 사용하게 된다. 

---

  오늘의 리눅스 커맨드는, 유사 bed format 의 라인들을 chromosome 별로 분리하기 위해 perl 코드를 작성하다가 구글링 해서 찾은 awk 커맨드 이다 .! 두둥 ! 아직 엄청 익숙하지는 않고 몇개 조건에 맞는 줄 출력하고 이럴 때 유용하게 쓰는데, 여러 파일로 분리해서 저장할 수도 있는 지점이 신기하고, 잘 알아두면 다른 경우에도 아주 유용할 것 같다.

awk '{close(f);f=$1}{print > f".bed"}'  [input bed]

  위 커맨드는 input bed 파일을 받아 첫 번째 컬럼의 이름으로 bed 파일을 만들고 읽어들인 라인들을 첫 번째 컬럼을 이름으로 가지는 bed파일로 direction 하여 출력하는... 내부적인 작동 방식은 모르겠고 저 커맨드를 치면 위의 유사 bed format 의 파일의 첫 번째 컬럼 (chromosome)별로 분리된 파일을 쇽쇽 생성해준다, 

  내가 사용한 데이터로 결과물을 만들면 아래와 같이 chromosome 별로 분리된 bed 파일을 예쁘게 얻을 수 있다. ".bed" 위치를 수정해서 원하는 확장자/파일이름 으로 만들어 줄 수도 있겠다.! 

 

awk 커맨드 실행해서 얻은 분리된 파일들 

 

Reference

seqanswers.com/forums/showthread.php?t=8115