[microbiome] long-read seq(nanopore, MinION) 16S rRNA analysis (진행중)

자내가 사용한 데이터는 microbiome 분석을 목적으로한 16s rRNA에 대한 MinION 시퀀스데이터이다. 

quality 에 통과한 fastq.gz를 barcode별로 받았고 taxanomic classification을 하는 것이 목적이다.

참 단순하게도 tool을 이용하는 것이지만 비주류 데이터를 쓰다보니 이 마저도 어려웠다.

그나마 다행인점은 아래 사용하는 모든 툴은 pip 이나 conda를 사용해 쉽게 설치할 수 있었다.

 

아래 두 개의 논문을 참고하여 파이프라인을 설계했다.

1. computational methods for 16s metabarcoding studies using nanopore sequencing data(link)

2. Benchmarking the MinION: Evaluating long reads for microbial profiling(link)

- 논문에 따르면 centrifuge> kraken2> kraken 순으로 정확도가 높다고 이야기하고 있으나, 그 차이가 크지 않다고 보인다.

- tabl1 을 확인하면 MinION의 quality 평균 12를 웃도는 것을 확인할 수 있었다. 이는 short read(평균30-40)와 다른 특징이다.

---

0. 환경 설정

가상환경 설정

-python 요구사항

centrifuge -> python[version='2.7.*|3.5.*|3.6.*|<3|>=3.6,<3.7.0a0|>=2.7,<2.8.0a0|>=3.5,<3.6.0a0|3.4.*']

 

conda create -n metagenome python=3.5

conda activate metagenome

 

필요한 프로그램

• nanofilt
• porechop
• centrifuge / kraken / epi2melabs

필요한 데이터

• bacterial 16S rRNA fasta file (NCBI refseq) : ftp://ftp.ncbi.nlm.nih.gov/refseq/TargetedLoci/Bacteria/bacteria.16SrRNA.fna.gz
• taxonomy file : https://ftp.ncbi.nlm.nih.gov/pub/taxonomy/taxdump.tar.gz
  • 여러 파일이 존재하는데 자세한 내용 참고 : README
  • nodes.dump 파일로부터 해당 Database에 몇 개의 species가 있는지 확인할 수 있다.
  • 2022년 NCBI refseq 기준, genus 105 242, species 2 023 060
  • 2024년 NCBI refseq 기준, genus 109 264 , species 2 118 227

 cat nodes.dmp | grep -w "genus" | cut -d "|" -f 1 | wc -l
 cat nodes.dmp | grep -w "species" | cut -d "|" -f 1 | wc -l

• seqence id to taxon id mapping file :kraken2-build --download-taxonomy --db $DBNAME
  • detail : 
  • cmd :
  • input
  • output:
• kraken DB index file : 3개의 .k2d 파일 

 

 

 

kraken2-build \

--add-to-library bacteria.16SrRNA.fna \

--db $DBNAME kraken2 \

-build --build --db $DBNAME

• centrifuge DB index file : 4개의 .cf 파일

centrifuge-build \

--conversion-table seqid2taxid.map \

--taxonomy-tree taxonomy/nodes.dmp \

--name-table taxonomy/names.dmp \

bacteria.16SrRNA.fna abv

---

1. removal low quality reads - NanoFilt

깃허브에 따르면 quality가 10 혹은 12이하 이거나 length가 500 이하이면 제거하는 과정이다. 이 기준이 적절한지는 좀 더 확인이 필요하다. (nanofilt의 대안이 되는 filtlong tool같은 경우 length가 1k 인 read는 모두 제거하고있다.)

+ 추가사항)  현재 진행하는 연구의 경우 품질보다 전체적인 조성이 중요하기 때문에, 이 단계는 건너 뛰고 진행했다.

 

설치 : conda install -c bioconda nanofilt

참고 자료 : 

[1] https://github.com/wdecoster/nanofilt

 

2. trimming of adaptor - porechop

conda install porechop

porechop은 2018년도에 업데이트를 중지했다.

하지만 유일무이한 nanopore 데이터에 대한 어댑터 제거하는 툴로 여전히 일반적으로 사용되는 것 같다.

illumina 데이터의 경우 Cutadapt 와 Trimmomatic 가 대표적으로 사용된다. 

default 파라매터를 실행하는 방법이다.

install : conda install -c bioconda porechop
cmd :porechop -i {FILE_DIR} -o {FILENAME}.fastq.gz | tee trimmed.log

참고 자료 : 

[1] porechop github : https://github.com/rrwick/Porechop

 

여기까지 진행하고 바로 centrifuge 를 돌렸는데, EPI2ME와 비교했을대 확연한 차이를 보였다. 문제가 있다.. minibar를 안돌렸기 때문인가.. -> centrifuge-kreport 명령어로 원하는 결과를 리포팅할 수 있었다. 또한 QC를 진행하지 않았기 때문이다.

 

3. sequence orientation - Minibar

nanopore 데이터의 경우 single-end 형식의 데이터이며, 내부적으로 foward, reverse strand 가 섞여있는 것이 특징이다.

방향(orientation)이 섞여있어 이를 맞춰주기 위해 이 같은 단계가 필요하다고 한다.

 

참고 자료 : 

[1] minibar git hub : https://github.com/calacademy-research/minibar

실제로 해보니 2번을 수동화 시킨것과 같다. adapter와 primer 정보를 인풋으로 넣어주면 이들을 제거하는 tool이다.

 

4. taxonomic assignmnet - centrifuge / kraken2

결과물이 안 좋다 뭐가 문제일까.. -> 데이터 베이스를 통일해주기로 했다.(ncbi 16s rRNA, bacteria)

비교시 필요한 기준 : # of total reads, # of classified reads, # of genus, genus rank

( # number 보다는 비율이 비교하기쉽다.)

EPI2ME 보여주는 결과를 정답이라고 생각했을 때,

kraken2 와 centrifuge 를 직접 돌린 결과와 편차가 많이 나는 것을 확인할 수 있었다. 

원인으로 1) 버전이 다른 ncbi refseq DB 데이터를 사용해서. 2) quality control 를 하지 않아서 라고 생각한다.

 

4.1 centrifuge

install : conda install -c bioconda centrifuge
cmd : centrifuge -x p_compressed+h+v -U {FASTQ_FILE} --report-file {FILENAME}.report.txt -S {FILENAME}.results.txt

에러 :

1. --verbose 옵션을 통해 어느 단계에서 에러가 났는지 확인한다.

2. 나 같은 경우 에러가 발생하지 않았지만, 빈 결과 파일이 만들어지는 문제가 있었다.

-x db_index_name 이 옵션을 절대 경로로 넣어줬기 때문이었다.

- command line을 수정하고 : /home/dayeon20/metaGenome/TOOLS/db_index/abv -> abv

- 현재 작업 위치에 abv.1.cf ~ abv.3.cf 파일을 옮겼다.

 

참고 자료 : 

[1] centrifuge manual(github) : https://github.com/DaehwanKimLab/centrifuge/blob/master/MANUAL\

[2] centrifuge-kreport code : https://github.com/jsh58/metagen/blob/master/centrifuge-kreport

[3] centrifuge home : http://ccb.jhu.edu/software/centrifuge/manual.shtml

https://github.com/DerrickWood/kraken2/wiki/Manual#inspecting-a-kraken-2-databases-contents

4.2 Kraken

레포트를 보면 clade 라는 개념이 나오는데 공통 조상에서 파생된 모든 자손을 의미한다.

install :conda install -c bioconda kraken2
cmd : kraken2 -db 16S_targetedLoci/ .{FASTQ_FILE} --output {FASTQ_FILE}.result.txt --report {FASTQ_FILE}.report.txt

참고 자료 : 

[1] kraken2 git hub : https://github.com/DerrickWood/kraken2/blob/master/docs/MANUAL.markdown

[2] 초보자를 위한 kraken2 tutorial : https://bioinformaticsworkbook.org/dataAnalysis/Metagenomics/Kraken.html#gsc.tab=0

[3] kraken2 wiki : https://github.com/DerrickWood/kraken2/wiki/Manual

4.3

MinION, Oxford 회사에서 자체적으로 tutorial를 제공하고 있다.

설치 : pip install  epi2melabs
cmd : 

참고 자료 : 

[1] post-EPI2ME 16S Analysis : https://labs.epi2me.io/notebooks/Analysis_of_EPI2ME_16S_CSV_Output.html

 

5. data visualization

https://genomics.sschmeier.com/ngs-taxonomic-investigation/index.html

https://github.com/marbl/Krona/wiki/Installing

 

6. calculate Diversity 

R.H. Whittaker 식물학자
the total species diversity in a landscape (gamma diversity) is determined by two different things, the mean species diversity in sites at a more local scale (alpha diversity) and the differentiation among those sites (beta diversity).

 

1. Alpha diversity

정의) 관심 지역의 평균 종 다양성
- 관심 지역의 크기는 상황에 따라 매우 다를 수 있으며, 특정 공간 규모에 제한되지 않음
- 생태계에서 생물 종들 간의 개체 수에 관계없이 얼마나 종의 수가 많은 가를 측정

- 종이 다양할 수록 (richness), 종들이 균등하게 분포할 수록(evenness) 종 다양성이 높다고 간주한다.

https://bioinfoblog.tistory.com/141

[Alpha diversity] Diversity metrics 비교 (Species richness, Species evenness)

공식)
- 생태학자들마다 조금씩 다른 정의의 alpha diversity를 사용

- Chao1, Shannon, Simpson, ACE, ...

- 미생물 분석에서 잘 사용하는 천랩 EZ biocloud 에서는
ACE, chao1, Simpson, NPShannon, Jackknife, Shannon 를 제공하고 있다.

 

1. Shannon

2. Shimson

 

- Beta diversity

- Gamma diversity

---

그 외, 각 단계별 여러 tool 존재한다.

https://github.com/B-UMMI/long-read-catalog/blob/master/README.md

---

+ 추가 자료)

- taxonomy ID 검색 : https://www.ncbi.nlm.nih.gov/Taxonomy/Browser/wwwtax.cgi?mode=info&id=1912

- 시퀀싱 분석의 모든 workflow : https://genomics.sschmeier.com/ngs-taxonomic-investigation/index.html

- phylogenetics 개념 : 

[1] https://ib.bioninja.com.au/standard-level/topic-5-evolution-and-biodi/54-cladistics/cladograms.html

 

[2] http://haasegen564s17.weebly.com/phylogenetics.html