微生物群集解析ツール QIIME 2 を使う．Part2.（QC，FeatureTable生成編）

QIIME 2 を使う．

QIIME2 はDNAのシーケンスデータから微生物解析を行うためのオープンソースのパイプラインである．クオリティの高いグラフや統計の処理を行うことが可能である．

今回は，前回の[微生物群集解析ツール QIIME 2 を使う．Part1.（インストール，ファイルインポート編）] (https://qiita.com/drafts/5d29bd8bbe1d19fdeec5/edit)の続きで，fastq ファイルのクオリティコントロール（QC）とFeature Table の作成を行う．

Fastq ファイルのクオリティコントロール（QC）を行う．

cutadapt にて、アダプター配列の除去。
341f (--p-front-f CCTAYGGGRBGCASCAG)と806r (--p-front-r GGACTACNNGGGTATCTAAT)を想定。

# cutadapt による341f (--p-front-f CCTAYGGGRBGCASCAG)と806r (--p-front-r GGACTACNNGGGTATCTAAT)のプライマー配列除去	
qiime cutadapt trim-paired --i-demultiplexed-sequences demux-paired-end.qza --p-cores 4 --p-front-f CCTAYGGGRBGCASCAG --p-front-r GGACTACNNGGGTATCTAAT --o-trimmed-sequences primer-trimmed-demux-paired-end.qza --verbose

# Summarise the reads
qiime demux summarize --i-data primer-trimmed-demux-paired-end.qza --o-visualization primer-trimmed-demux-paired-end.qzv

DADA2はIllumina のアンプリコン配列を検出し修正するパイプラインである．Phix配列やキメラ配列をフィルタリングしてくれる．

公式チュートリアルでは，dada2 denoise-singleになっているが，今回のデータはpaired-end のデータなのでdada2 denoise-pairedで実行する．
deblur error "Argument to parameter 'demux' is not a subtype of SampleData[SequencesWithQuality]."

dada2 denoise-paired メソッドでは，4つのパラメータを与える．

--p-trim-left-f x：各シーケンスのforward 配列の初めから塩基長xまでを削除する．
--p-trim-left-r y：各シーケンスのReverse 配列の初めから塩基長yまでを削除する．

出力されたクオリティのプロットを見て、クオリティが20以下の部分をカットしつつ。
採用する最大配列長には結合したペアエンドリードが十分重なる長さを指定することが求められる．

--p-trunc-len-f m：各シーケンスのforward 配列から採用する最大配列長mを決定する． mを超える塩基長は除去する．
--p-trunc-len-r n：各シーケンスのreverse 配列から採用する最大配列長nを決定する． nを超える塩基長は除去する．

v3-v4領域だと460bpなので、最低でも 20bp 以上のオーバーラップ領域を加味して、270,210とか480を超える長さで実行する．

--p-n-threads k：解析に使用するスレッド数を指定する．

プライマー配列は取り除く必要があるので、先頭側のカット長はプライマー配列長を通常指定する．

qiime dada2 denoise-paired \
  --i-demultiplexed-seqs primer-trimmed-demux-paired-end.qza \
  --p-trim-left-f 0 \
  --p-trim-left-r 0 \
  --p-trunc-len-f 280 \
  --p-trunc-len-r 220 \
  --o-representative-sequences rep-seqs-dada2.qza \
  --o-table table-dada2.qza \
  --o-denoising-stats stats-dada2.qza \
  --p-n-threads 30
# 名前を変更する（別にしなくても良いがqiime2のチュートリアルでの汎用性を考えてここでは変更しておく）
mv rep-seqs-dada2.qza rep-seqs.qza
mv table-dada2.qza table.qza

Feature Table の生成とサマライズを行う．

クオリティのフィルタリングが完了したら，結果を調べるためにFeature Table を生成する．

feature-table summarize:各サンプルに関連付けられているシーケンスの数、各特徴、それらの分布のヒストグラムを出力する．
feature-table tabulate-seqs:Feature ID のシーケンスへのマッピングとNCBI database に対するBLAST検索のリンクを提供する．

# metadata ファイルをダウンロードする
wget -O "sample-metadata.tsv" "https://data.qiime2.org/2017.12/tutorials/moving-pictures/sample_metadata.tsv"
#or
curl -sL "https://data.qiime2.org/2017.12/tutorials/moving-pictures/sample_metadata.tsv" > "sample-metadata.tsv"

# サマライズの実行
qiime feature-table summarize \
  --i-table table.qza \
  --o-visualization table.qzv \
  --m-sample-metadata-file sample-metadata.tsv

qiime feature-table tabulate-seqs \
  --i-data rep-seqs.qza \
  --o-visualization rep-seqs.qzv

補足

DAD2の公式ドキュメントにおいて、「How can I remove primers/adapters/etc from my amplicon reads?」という項目にて、重要なことが書かれている。

プライマー配列の除去について、より複雑な状況（例えば、不均一性スペーサー、ITS領域などの可変長アンプリコン）の場合は、cutadaptツールまたはtrimmomaticツールをお勧めします。
プライマー配列が残っていると後の解析に影響してしまうので、しっかり除去しよう。

If primers are at the start of your reads and are a constant length (a common scenario) you can use the trimLeft = c(FWD_PRIMER_LEN, REV_PRIMER_LEN) argument of dada2’s filtering functions to remove the primers.

For more complex situations (e.g. heterogenity spacers, variable length amplicons such as the ITS region) we recommend the cutadapt tool or the trimmomatic tool. You can see a worked example, including verification of primer removal, in the ITS-specific version of the DADA2 workflow.

Please double-check that your primers have been removed! The ambiguous nucleotides in unremoved primers will interfere with the dada2 pipeline.