Entrez API를 사용하려면 NCBI Insights에서 제공하는 문서인 How to Download Bacterial Genomes Using the Entrez API를 참조하기 바란다.
여기서 소개하는 방법은 이제 더 이상 권장하지 않는다. 왜냐하면 ncbi-genome-download라는 탁월한 스크립트가 있기 때문이다. 기관의 내부 전산망을 이용하는 경우 SSL 인증서 문제로 약간 애를 먹을 수는 있을 것이다. — Haeyoung Jeong 2021/06/09 14:37
[최신 소식] 2016년 9월 22일을 기하여 ftp://ftp.ncbi.nlm.nih.gov/genomes/all 하위 구조가 바뀌었다. 한 디렉토리 안에 너무나 많은 서브디렉토리가 존재하는 문제를 피하기 위하여 assembly accession prefix인 9자리 숫자를 3개씩 끊어서 3단계의 디렉토리 구조로 전환한 것이다(원본 문서).
the data currently in genomes/all/GCA_000001405.23_GRCh38.p8
will be moved to genomes/all/GCA/000/001/405/GCA_000001405.23_GRCh38.p8
the data currently in genomes/all/GCF_001696305.1_UCN72.1
will be moved to genomes/all/GCF/001/696/305/GCF_001696305.1_UCN72.1
그러나 ftp://ftp.ncbi.nlm.nih.gov/genomes/(archive|genbank|refseq)는 아직 현재의 체제를 유지하고 있다. Genome download(FTP)와 관련한 일반적인 질문과 답변은 FAQ 사이트를 참조하라. 이 FAQ 사이트에는 대단히 중요한 정보가 많으니 일독을 권한다.
Entrez 창에 다음과 같은 질의어를 넣어서 나타나는 결과물 중 <Genomes> 섹션의 Assembly에서 원하는 것을 찾는다. 종합적인 정보를 알고 싶다면 Genome을 클릭한다. 이 문서는 Entrez 질의어 구성법에 관한 것이 아니니 자세한 설명은 하지 않는다.
"Escherchia coli"[organism] AND "BL21(DE3)"[infraspecific name]
위 질의어에서 BL21(DE3)를 따옴표로 둘러싼 첫은 괄호를 특별한 기호로 인식하기 때문이다. 결과가 궁금하면 여기를 클릭해 보라. Assembly 검색 결과는 총 4건으로 나오지만 실제로는 2건이다. 각각에 대하여 GenBank와 RefSeq 레코드가 별도로 존재하므로(유전체 서열은 똑같다) 4건으로 나오는 것이다. 사실 이렇게까지 잘 짜 넣지 않고 대충 콩떡같이 질의어를 던져도 Entrez는 찰떡같이 해석하여 결과를 잘 찾아낸다.
Escherichia coli "BL21(DE3)"
이렇게 물어본 것은 실제로 내부적으로 다음과 같이 처리한다. 모든 레코드가 Organism: AAAA BBBB, Infraspecific name: Strain: CCCC의 표준 규격을 잘 맞추어서 제출되었을 것 같지는 않다. 그래서 아무리 느슨하게 질의어를 넣어도 최대한의 검색 결과를 끌어내려는 메커니즘이 돌아가는 것으로 보인다. latest.. 이후는 어떤 검색어를 넣어도 자동적으로 붙는 것 같다.
("Escherichia coli"[Organism] OR Escherichia coli[All Fields]) AND "BL21(DE3)"[All Fields] AND (latest[filter] AND all[filter] NOT anomalous[filter])
만약 species의 NCBI taxonomy ID를 안다면 웹브라우저 주소창에 다음과 같이 URL을 입력하여 genome sequencing 현황을 볼 수 있다. 167은 E. coli의 taxonomy ID라서 상당히 많은 결과가 나올 것이다.
https://www.ncbi.nlm.nih.gov/genome/genomes/167
아무리 웹 화면이 편리하다 한들 원하는 모든 정보에 대한 답을 담고 있지는 않다. '오늘 날짜 기준으로 RefSeq에는 총 몇 건의 유전체 정보가 존재하고, 이 중에서 Microcystis aeruginosa의 완성된 유전체 어셈블리만을 선별하여 마우스로 클릭하지 않고 한번에 전체를 다운로드할 수는 없을까?' 물론 가능하다. 그러려면 ftp site에서 관련 정보 파일(= assembly_summary.txt)을 먼저 받아서 조건에 맞는 라인을 찾은 뒤, 다운로드 주소를 일괄적으로 확보하여 wget 또는 curl로 받으면 된다. 라인 하나는 하나의 assembly에 대한 정보를 수록한다. assembly_summary.txt는 다음과 같이 주요 kingdom(archaea, bacteria, fungi 등) 별로 따로 존재한다. GenBank 카테고리의 파일을 쓰려면 아래의 주소에서 refseq을 전부 genbank로 치환하면 된다. RefSeq 자료는 NCBI의 정책에 따라 일관된 방식으로 genome annotation을 실시해 놓는다. 그러나 GenBak의 자료가 항상 100% RefSeq에 있는 것은 아니니 목적에 맞게 적당히 선택하여 사용하도록 한다. 2018년 2월 12일 현재 확인한 바에 의하면 Paenibacillus의 경우 RefSeq의 자료가 20건 적었다. 이것은 최근에 제출된 자료가 아직 RefSeq로 넘어오지 않은 데서 기인한 것이 아니다! 왜 그런지는 추후에 자세히 알아보겠다.
Klebsiella pbeumonia 중에서 Complete Genome 혹은 Chromosome 수준의 것만을 모으는 예제(genome nucleotide sequence file)를 가지고 설명한다. 각 과정은 다음과 같다.
받은 파일의 이름은 (assembly_accession)_(asm_name)_(나머지 부분)의 형식을 갖는다. 예를 들자면 GCF_002264385.1_ASM226438v1_genomic.fna와 같은 방식이다. 균주의 이름이나 우리가 논문 등에서 흔히 접하는 서열 accession number(CP######, NC_######)는 아예 나타나지도 않고 asm_name(ASM#####v#)은 unique하지도 않다. 따라서 파일의 이름만 보고도 대략적인 정보를 파악할 수 이도록 이를 적절히 수정해야 한다. 새로운 파일명에 사용할 모든 정보는 당연히 list_all에 있으므로 이를 적절히 조합하여 만든다. 이 과정에 도움을 주는 스크립트가 custom Perl script인 my_rename1.pl과 my_rename2.pl이다(링크). 혹시 파일명에 부적합한 값(예: 공백)이 있을지도 모르니 눈으로 한번 확인하는 것이 좋다. 이를 다운로드한 파일과 연결할 키는 (assembly_accession)_(asm_name)이다(my_rename1.pl의 결과 파일인 id2name의 첫번째 컬럼). (assembly_accession)만을 키로 사용해도 될 것이다.
$ curl -o refseq_bacteria_assembly_summary_`date +%Y-%m-%d`.txt ftp://ftp.ncbi.nlm.nih.gov/genomes/refseq/bacteria/assembly_summary.txt
# 정규식 안의 ^는 없어도 괜찮다.
$ awk -F "\t" '$8~/^Klebsiella pneumoniae/ && $12~/^Complete Genome|^Chromosome/' refseq_bacteria_assembly_summary_2017-07-04.txt > list_all
$ awk -F "\t" '{print $20}' list_all > list_all_ftp_path
# fna, gbff, or faa file? 용도에 맞게 커맨드 라인을 작성할 것.
$ sed -r 's|(ftp://ftp.ncbi.nlm.nih.gov/genomes/all/.+\/)(GCF_.+)|\1\2\/\2_genomic.fna.gz|' list_all_ftp_path > list_all_files_with_path_fna
$ sed -r 's|(ftp://ftp.ncbi.nlm.nih.gov/genomes/all/.+\/)(GCF_.+)|\1\2\/\2_genomic.gbff.gz|' list_all_ftp_path > list_all_files_with_path_gbff
$ sed -r 's|(ftp://ftp.ncbi.nlm.nih.gov/genomes/all/.+\/)(GCF_.+)|\1\2\/\2_protein.faa.gz|' list_all_ftp_path > list_all_files_with_path_faa
# ~_protein.faa.gz와 ~_translated_cds.faa.gz는 서열 ID와 description의 모양새가 다름에 유의하라. locus tag등 월등하게 많은 정보를 수록한 것은 후자이다.
$ mkdir 01_fna
$ cd 01_fna
# GNU parallel을 응용하여 여러 파일을 동시에 받을 수도 있을 것이다.
$ cat ../list_all_files_with_path_fna | while read f; do wget $f; done
# 또는 cat ../list_all_files_with_path_fna | parallel -j6 wget
$ gzip -d *gz
$ my_rename1.pl ../list_all > id2name
# id2name 파일의 두번째 컬럼은 파일 이름으로 쓰일 것이므로 특수문자가 있거나 중복되면 안된다.
# 공백은 밑줄('_')로 치환된다.
# 수정을 하였다면 id2name.mod로 변경한다.
$ ls *fna | while read f; do my_rename2.pl id2name $f; done
$ rm -f GCF*
id2name 파일을 매개로 하여 GCF_000469945.2_ASM46994v2_genomic.fna.gz라는 파일이 Paenibacillus_sp._P22_GCF_000469945.2.fna로 바뀌게 된다. 변경 후 파일명의 앞부분을 차지하는 균주 이름은 당연히 assembly summary file의 여덟번째 컬럼에서 유래한다. 하지만 실제 assembly summary 파일을 열어보면 파일명으로 쓰기에 적합하지 않는 문자가 들어있는 경우가 상당히 많다. 예를 들어 공백, 세미콜론, 괄호 등이 그러하다. 이를 되도록 철저하게 밑줄로 바꾸는 것이 중요하다. 이미 파일 이름을 다 바꾼 상태에서 GCF_000469945.2.fna처럼 assembly accession만 남기고 싶다면 [정해영의 블로그] Bash, shell: 구분자가 포함된 문자열의 조작 기법을 참조하라.
$ xargs wget < ../list_all_files_with_path_fna
multiple argument를 다루려면 어떻게 해야 하는가? 조금 더 알아봐라. 인수 위치의 지정은 다음과 같이 한다. 하지만 이 사례에서는 하나의 gff file에 모든 출력이 다 기록된다는 문제가 있다.
$ ls *fna | xargs -I{} barrnap --threads 16 {} > {}.rna.gff