Microbiome Database를 만들어볼까? -NCBI편 3-

지난 시간에 이곳에서 귀한 자료를 받아 봤을 겁니다.
그게 무엇이냐!!!

지금까지 공개된 bacteria의 서열들이죠 정확하게 말하면
현재까지 수많은 연구자들이 자발적(자의타의)으로
공개해준 것을 NCBI가 아름답게 정리한 RNA서열들만 다운받았습니다.

근데 우리는 RNA 서열이 아니라 rRNA 서열이 필요하죠

이전 글에서 우리는 rna라는 폴더 안에 12만개에 달하는 파일들을 다운받았습니다.
아마 그냥 ls하시면 어쩌구 저쩌구 long 할겁니다.
ls로 불러오기에 item이 너무 많다 이거죠
그럼 어쩌지?? @.@

우리에겐 귀도 훃님의 파이썬이 있지 않겠습니꽈

간단하게 다음과 같은 스크립트를 뚝딱 뚝딱 만들어보죠

import os,glob
from Bio import SeqIO
for files in glob.glob('rna/*gz`):
    for rec in SeqIO.parse(gzip.open(files), format='fasta'):
        name =rec.description
        seq = rec.seq
        if name.find('[product=16S ribosomal RNA]') != -1:
            print '>%s\n%s\n'.format(name,seq)

(python2.7에 Biopython이 설치되어 있어야 하고 *_rna_form_genomic.fna.gz 파일이 rna폴더 밑에 위치하고 있어야 합니다.)

$python script.py > ncbi.fa

라고 해주면 헤더에 "[product=16S ribosomal RNA]"가 포함된 aka 우리가 원하는 바로 그것! 16S rRNA 서열을 각 종에서 샤샤샥 ncbi.fa라는 파일에 저장할 수가 있습니다.

-12만개 파일 읽어 오는거라 순차적으로 하면 3-4시간 걸릴것이고
subprocess로 잘 해주시면 적어도 십수분? 30분이내면 충분히 끝날각 되겠습니다.

출처: SM Town