페이지네임을 한글로 길게 만들었더니 저장이 되지 않는 현상을 발견하였다. DoluWiki의 문제인 듯.
NCBI에서 제공하는 공식 문서 Genome Download (FTP) FAQ에 유용한 정보가 많으므로 이것을 일단 숙독하는 것이 좋을 것이다. 특히 'How can I find the sequence and annotatio of my genome of interest?' 질문이 중요하다. 동일한 유전체 서열에 대하여 왜 RefSeq과 GenBank의 레코드가 따로 존재하는지, 그리고 일부 유전체는 왜 GenBank에만 존재하는지 등에 대해서도 그 이유를 알아야 한다.
간단한 방법으로는 ncbi-genome-download가 편리하다. 개별적인 assembly accession은 물론 사용자가 정의한 유전체의 그룹 단위로 다운로드할 수 있어 매우 유용하다. 맨 마지막 옵션으로는 bacteria, archaea, fungi, viral과 같은 group 정보를 반드시 주어야 한다. 이 스크립트는 실행할 때 날짜 기준으로 $HOME/.cache/ncbi-genome-download에 assembly summary file의 캐시를 생성한다. –no-cache 옵션을 사용하면 이를 무시할 수 있다. '–assembly-accession'('-A'와 동일) 파라미터는 단일 accession, 쉼표로 구분한 accession, 또는 파일 형태로 수록한 accession을 전부 받아들일 수 있어 매우 유용하다. 단지 종 및 균주 명칭만을 입수하기 위함이라면 뒤에서 소개할 NCBI EDirect 명령어 조합(esearch | efetch | xtract)보다 훨씬 빠르다.
$ ncbi-genome-download bacteria # RefSeq의 모든 박테리아 유전체 다운로드
$ ncbi-genome-download --section genbank fungi # GenBank에서 모든 진균 유전체 다운로드
# 실제 다운로드는 하지 않고 accession number와 간단한 정보만 출력('--dry-run' or '-n')
$ ncbi-genome-download --dry-run --genus "Streptomyces" bacteria
$ ncbi-genome-download --type-materials type,reference
$ ncbi-genome-download --format fasta --taxid 511145 bacteria
$ ncbi-genome-download --assembly-accessions GCF_000146875.3,GCF_000833145.1 bacteria # 소수점 이하(버전)까지 지정해야 함
Accession number를 이용하여 뉴클레오티드 혹은 단백질 서열을 받으려면 ncbi-acc-download. 다운로드하려는 분자 유형이 단백질이라면 -m 또는 -molecule {nucleotide,protein} 옵션으로 지정해야 된다. 그러나 이 방법으로는 assembly accession(예: GCA_011389495.1 또는 GCF_011389495.1)1)을 사용한 다운로드는 불가하다. 다음 장에서 설명할 NCBI E-Direct 유틸리티도 PubMed나 서열 또는 유전체 정보의 검색 및 일괄 다운로드에 편리하게 사용할 수 있다.
위에서 소개한 두 종류의 유틸리티는 분류학적 정보나 assembly accession number를 근거로 하여 유전체 정보를 일괄적으로 다운로드하는 데에 적합하다. 그러나 NCBI에서는 등록된 모든 유전체 염기서열에 대하여 균주의 명칭, accession number, 제출자, 등록일, 균주 출처 등 스무 가지가 넘는 부가 정보를 체계적으로 정리하여 다음의 위치에서 텍스트 파일(assembly_summary_DB.txt)로 제공하므로, 여기에서 다운로드할 수 있는 파일의 FTP path를 찾아내어 wget 명령으로 일괄 다운로드하는 것이 훨씬 편리하다. 다음 두 가지의 assembly summary file은 모든 등록된 모든 생명체의 유전체 조립 정보를 GenBank 및 RefSeq로 나누어서 수록한 것이다. assembly file 자체와 각 컬럼의 의미에 대한 설명이 필요하면 README 파일을 참조한다.괄호 안의 숫자는 각각 2021년 8월 5일과 2023년 6월 22일 기준의 파일 크기이다. 약 2년 사이에 유전체 DB의 엔트리가 얼마나 증가했는지 알 수 있다.
RefSeq에 등록된 박테리아의 유전체 조립 정보를 참고하고 싶다면 다음의 파일을 활용한다. RefSeq가 아닌 GenBank 자료를 원한다면 URL의 'refseq'을 'genbank'로 치환하면 되고, 'bacteria'를 'archaea', 'fungi', 'invertebrate' 등으로 치환하면 다른 서브셋에 대한 정보를 가져오게 된다.
https://ftp.ncbi.nlm.nih.gov/genomes/refseq/bacteria/assembly_summary.txt
우선 오늘 날짜의 RefSeq assembly summary report 파일(bacteria)을 다운로드한다. Curl 대신 wget을 사용해도 된다. 매일 데이터가 갱신되므로 다운로드한 assembly summary report 파일의 이름에 날짜를 포함하는 것은 나름대로의 의미가 있다.
$ curl -o assembly_summary_refseq_bacteria_`date +%Y-%m-%d`.txt \ https://ftp.ncbi.nlm.nih.gov/genomes/refseq/bacteria/assembly_summary.txt
특정 종(예: Klebsiella pneumoniae)에 대하여 assembly level = 'Complete genome' 혹은 'Chromosome'으로 구성된 유전체 정보에 해당하는 라인만을 선택하여 FTP path 정보(20번째 컬럼)를 추출해 본다. Assembly level을 기준으로 선별하지 않으면 데이터의 분량이 크게 늘어날 것이다.
$ awk -F "\t" '$8~/^Klebsiella pneumoniae/ && $12~/^Complete Genome|^Chromosome/' bacteria_assembly_summary_refseq_*.txt > list_all
$ awk -F "\t" '{print $20}' list_all > list_all_ftp_path
이어서 다운로드할 파일의 URL을 완성한다. sed 치환 뒤의 파일명에서 'genomic.fna.gz'은 실제 다운로드하고 싶은 파일에 맞추어 수정한다(예: genomic.fna, cds_from_genomic.fna, genomic.gbff, protein.faa 등). 만약 list_all 파일을 GenBank에서 받은 assembly summary 파일에서 만들었다면, 다음 명령어에서 GCA를 GCF로 바꾸면 된다.
$ 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
새 디렉토리를 만들어 이동한 뒤 wget과 목록 파일(list_all_files_with_path_fna)을 이용하여 파일을 다운로드한다.
$ mkdir 01_fna; cd 01_fna $ cat ../list_all_files_with_path_fna | while read f; do wget $f; done $ gzip –d *gz
GNU parallel이 설치되어 있다면 한번에 지정한 수(-j6)의 파일을 동시에 다운로드할 수 있다.
$ cat list_all_files_with_path_fna | parallel -j6 wget
Parallel을 사용하면 간혹 접속 상태가 나쁠 때 다음 파일의 다운로드로 건너 뛰기도 한다. 이를 점검하려면 다음과 같이 표준입력과 에러를 구조화하여 출력하도록 parallel을 실행한 다음 표준에러 파일을 한데 모아서 문서편집기로 확인해야 한다. 'saved'로 표시된 파일의 수가 list_all_files_with_path_fna의 라인 수와 같아야 한다. 'failed' 메시지가 나온 것이 무엇인지 알아내어 다시 다운로드를 시도한다.
$ cat list_all_files_with_path_fna | parallel -j6 --results outdir wget
$ find outdir/ -name stderr -exec cat {} \; > stderr.all
$ grep -B 4 failed stderr.all
다운로드한 파일의 이름은 다음과 같은 형태라서 균주 정보를 직접적으로 유추하기가 어렵고 너무 길다.
GCF_000001405.38_GRCh38.p12_genomic.fna.gz GCF_000001215.4_Release_6_plus_ISO1_MT_genomic.fna.gz
따라서 나머지 필드는 다 제거하고 assembly accession(GCF_#########.#)만 남겨서 파일 이름이 되게끔 하는 것이 바람직할 것이다.
$ for i in $(ls *.fna.gz) > do > mv $i $(cut -d_ -f1,2 <<<$i).fna > done # GCF_000001405.38_GRCh38.p12_genomic.fna.gz -> GCF_000001405.38.fna.gz로 변경
또는 학명과 균주 이름으로 파일이름을 바꾸는 것도 좋다. 이를 위해서 my_rename1.pl과 my_rename2.pl를 활용해 보자. 먼저 my_rename1.pl을 사용하여 id2name 파일을 만든다.
$ my_rename1.pl ../list_all > id2name
my_rename1.pl 스크립트는 list_all 파일(assembly summary에서 선별된 row로 이루어진 부분집합)의 8번 컬럼(organism_name; NCBI taxid를 기준으로 가져온 값)과 9번 컬럼(infraspecific_name)을 우선적으로 처리하여 학명 + strain에 해당하는 정보를 우선적으로 생성하며, 적절한 species 하위 식별자가 없는 경우에는 16번 컬럼(asm_name)을 활용하도록 만들었다. 8번 컬럼에는 species 명칭만 있는 것이 현재의 원칙이지만, 과거에 등록된 정보 중에는 strain level까지 정의된 것이 있어서 스크립트 내에서 처리하기가 번거로운 것은 사실이다. NCBI BioSample 데이터베이스의 ‘Organism’ 정보 역시 assembly summary의 8번 컬럼과 비슷한 상황이다. 이를 활용하는 방법은 후술할 NCBI E-Direct 활용법을 참고하기 바란다.
id2name 파일의 두 번째 컬럼이 변경된 파일의 이름 base로 쓰일 것이므로 공백이나 특수문자 등 적합하지 않은 문자가 없도록 적절히 편집하여 id2name.mod로 저장한다. 괄호 표시도 밑줄로 대체하는 것이 바람직하다. 간혹 동일한 균주를 서로 다른 곳에서 2회 이상 시퀀싱하여 등록한 것을 발견하게 된다. 이런 것은 바뀐 후의 파일 이름이 동일해져 버리므로, assembly 상태가 좋은 것만을 남기는 것이 좋다. 다음을 실행하면 수정이 완료된 FASTA file을 얻는다.
$ ls *fna | while read f; do my_rename2.pl id2name.mod $f; done