Taxonomy annotation for amplicon sequencing - Biomedical Data Science Laboratory

Taxonomy annotation

소개

미생물의 유전자는 진화 과정을 통해 끊임 없이 변화한다. 동시에, 선조의 유전 정보 혹은 서열은 후대로 이어진다. 진화적으로 근연 관계에 있는 미생물은 비슷한 유전자 서열을 가지고 있으며, 역으로 유전자 서열을 이용하여 미생물의 계통을 예측할 수 있다. rRNA는 모든 미생물이 가진 유전자로, 이 유전자의 서열을 이용해 미생물의 taxon(분류군)을 유추할 수 있다.

이러한 예측을 위해서는 미생물의 (1) rRNA 서열과 (2) taxon 에 대한 데이터가 필요하다. 이 데이터를 이용하여 rRNA 서열로 부터 taxon을 예측하는 모델을 구성할 수 있고, 그 모델을 이용하여 서열 정보로 부터 taxon 정보를 유추할 수 있다. 여기서는 QIIME2 package를 이용한 taxon 예측 방법을 설명한다. 그리고, 이어서 서열로 부터 taxon을 예측하는 Naive Bayesian 모델을 학습하는 방법을 설명한다.

Annotate Taxonomy

Taxonomy를 예측할 수 있는 classification model을 이용하여 taxonomy를 예측할 수 있다. Classification model은 rRNA의 서열 – taxonomy 쌍 정보를 이용하여 구성할 수 있으며, Naive Bayesian model이 좋은 예측 성능을 보여 널리 활용된다. rRNA-taxonomy database로는 SILVA나 greengene 가 널리 이용된다. Classification model이 확보되면 아래 Qiime2 명령어를 이용하여 taxonomy를 예측한다.

qiime feature-classifier classify-sklearn \
        --i-classifier silva-138.1-ssu-nr99-v3v4-classifier.qza \
        --i-reads sequences.qza \
        --o-classification taxonomy.qza

Construct taxonomy classifier

QIIME2 package와 RESCRIPt plugin을 이용하여 Naive Bayesian classification model를 구성하는 방법을 설명하며, 특히, 16S rRNA V3-V4 region을 대상으로 하는 classifier 구축하는 방법을 설명한다.

Taxonomy prediction을 위해서 SILVA database를 이용한다. SILVA는 미생물의 rRNA 서열과 taxon 정보를 포함하는 데이터베이스이다 (SILVA). 오래동안 greengene database가 널리 이용되었으나, 최근 자료를 반영하지 못하는 단점이 있어 SILVA database가 널리 활용된다.

여기서 설명하는 내용은 아래의 링크의 내용에 바탕 하여 작성하였다. 이 내용은 Conda를 통해 QIIME2 package를 설치하였음을 전제로 한다.

Based on: https://forum.qiime2.org/t/processing-filtering-and-evaluating-the-silva-database-and-other-reference-sequence-data-with-rescript/15494

Step 1. Install RESCRIPt within QIIME2 environment

SILVA database로 부터 데이터를 받고, 데이터를 정리하고 위해서 RESCRIPt를 이용한다. RESCRIPt는 SILVA와 같은 sequence database를 관리하기 위한 Python package이다 (RESCRIPt). 여기서는 QIIME2 package내에 RESCRIPt plugin을 설치하는 과정을 설명한다.

아래 명령어 실행을 통해 QIIME2 환경을 활성화 시키고, 필요한 library를 설치하고, RESCRIPt package를 설치한다.

conda activate qiime2-2022.8

# Installation of libraries required for RESCRIPt
conda install -c conda-forge -c bioconda -c qiime2 -c defaults xmltodict

# Installation of RESCRIPt
pip install git+https://github.com/bokulich-lab/RESCRIPt.git

Step 2. Download SILVA database

SILVA database로 부터 데이터를 받기 위해 아래 명령어를 수행한다. Version 138.1을 받으며 rRNA small subunit (SSU)에 대한 서열을 받는 명령어이며, 99% sequence identify로 서열을 clustering한 결과이다.

qiime rescript get-silva-data \
    --p-version '138.1' \
    --p-target 'SSURef_NR99' \
    --p-include-species-labels \
    --o-silva-sequences silva-138.1-ssu-nr99-rna-seqs.qza \
    --o-silva-taxonomy silva-138.1-ssu-nr99-tax.qza

Step 3. Preprocess SILVA database

다운 받은 database의 데이터의 질을 점검하고, 낮은 품질의 데이터를 삭제한다.

Convert RNA sequences into DNA sequences
Screen out low-quality sequences
Filter out short sequences
Remove redundant sequence-taxonomy pairs

# 1. Convert RNA sequence into DNA sequence
qiime rescript reverse-transcribe \
    --i-rna-sequences silva-138.1-ssu-nr99-rna-seqs.qza \
    --o-dna-sequences silva-138.1-ssu-nr99-seqs.qza

# 2. Culling low-quality sequences
qiime rescript cull-seqs \
    --i-sequences silva-138.1-ssu-nr99-seqs.qza \
    --o-clean-sequences silva-138.1-ssu-nr99-seqs-cleaned.qza \
    --p-n-jobs 20

# 3. Filter sequence by minimum length
qiime rescript filter-seqs-length-by-taxon \
    --i-sequences silva-138.1-ssu-nr99-seqs-cleaned.qza \
    --i-taxonomy silva-138.1-ssu-nr99-tax.qza \
    --p-labels Archaea Bacteria Eukaryota \
    --p-min-lens 900 1200 1400 \
    --o-filtered-seqs silva-138.1-ssu-nr99-seqs-filt.qza \
    --o-discarded-seqs silva-138.1-ssu-nr99-seqs-discard.qza

# 4. Remove redundant sequence
qiime rescript dereplicate \
    --i-sequences silva-138.1-ssu-nr99-seqs-filt.qza  \
    --i-taxa silva-138.1-ssu-nr99-tax.qza \
    --p-mode 'uniq' \
    --o-dereplicated-sequences silva-138.1-ssu-nr99-seqs-derep-uniq.qza \
    --o-dereplicated-taxa silva-138.1-ssu-nr99-tax-derep-uniq.qza \
    --p-threads 20

Step 4. Extract V3-V4 region

여기서는 V3-V4 region을 대상으로 amplicon sequencing한 경우를 가정한다. 따라서 classification model을 구성할 때도 해당 부위의 서열 만을 이용한다.

Extract V3-V4 region
Remove redundant sequence-taxonomy pairs

# 1. Extract V3-V4 region
qiime feature-classifier extract-reads \
    --i-sequences silva-138.1-ssu-nr99-seqs-derep-uniq.qza \
    --p-f-primer CCTACGGGNGGCWGCAG \
    --p-r-primer GACTACHVGGGTATCTAATCC \
    --p-n-jobs 10 \
    --p-read-orientation 'forward' \
    --o-reads silva-138.1-ssu-nr99-seqs-v3v4.qza

qiime rescript dereplicate \
    --i-sequences silva-138.1-ssu-nr99-seqs-v3v4.qza \
    --i-taxa silva-138.1-ssu-nr99-tax-derep-uniq.qza \
    --p-mode 'uniq' \
    --o-dereplicated-sequences silva-138.1-ssu-nr99-seqs-v3v4-uniq.qza \
    --o-dereplicated-taxa  silva-138.1-ssu-nr99-tax-v3v4-derep-uniq.qza

Step 5. Build classifier

마지막으로 classificatio model을 구성한다. QIIME2 의 Naive Bayesian classifier 모델을 이용하며, SILVE nr99의 V3-V4 부위 데이터를 이용한다.

qiime feature-classifier fit-classifier-naive-bayes \
    --i-reference-reads silva-138.1-ssu-nr99-seqs-v3v4-uniq.qza \
    --i-reference-taxonomy silva-138.1-ssu-nr99-tax-v3v4-derep-uniq.qza \
    --o-classifier silva-138.1-ssu-nr99-v3v4-classifier.qza