### 方法1. DADA2   feature table   双端合并，去除嵌合体，截去接头序列降噪生成feature table 一步完成
cd $workdir
mkdir dada2

echo dada2 start
date
cd dada2
	
	qiime dada2 denoise-paired \
	  --i-demultiplexed-seqs $workdir/demux.qza \
	  --p-n-threads 1 \
	  --p-trim-left-f 29 --p-trim-left-r 18 \
	  --p-trunc-len-f 0 --p-trunc-len-r 0 \
	  --o-table dada2-table_biom.qza \
	  --o-representative-sequences dada2-repset-seqs.qza \
	  --o-denoising-stats denoising-stats.qza

echo dada2 end
date

### 方法2. 外部导入特征表和代表序列(常用)  
cd $workdir
mkdir import_table
cd import_table
# 上传我们生成的OTU表otu_table.txt和代表序列rep_set.fa
# 转换文本为Biom1.0，注意biom --version 2.1.5/8可以，2.1.7报错
qiime tools import --input-path $workdir/otu_table.biom \
  --type 'FeatureTable[Frequency]' --input-format BIOMV100Format \
  --output-path table_biom.qza

qiime tools import --input-path $workdir/otus.fa \
  --type 'FeatureData[Sequence]' \
 --output-path repset-seqs.qza

### 方法3. deblur  feature table    deblur只能输入合并好的序列，  速度比dada2快
#因此先做序列处理，去除引物，合并，过滤低质量等等
# Trim amplicon primers  去除引物序列   # 341f   # 806r
qiime cutadapt trim-paired \
  --i-demultiplexed-sequences $workdir/1.import_data/demux.qza \
  --p-cores 4 \
  --p-no-indels \
  --p-front-f CCTAYGGGRBGCASCAG \
  --p-front-r GGACTACNNGGGTATCTAAT  \
  --o-trimmed-sequences primer-trimmed-demux.qza 

# Summarise the reads  查看结果
qiime demux summarize \
  --i-data primer-trimmed-demux.qza \
  --o-visualization primer-trimmed-demux.qzv \


#双端合并 

qiime vsearch join-pairs \
  --i-demultiplexed-seqs primer-trimmed-demux.qza \
  --p-threads  4 \
  --o-joined-sequences demux-joined.qza
#查看合并结果
qiime demux summarize \
  --i-data demux-joined.qza \
  --o-visualization demux-joined.qzv

#质控过滤数据  
qiime quality-filter q-score \
  --i-demux demux-joined.qza \
  --o-filtered-sequences demux-joined-filtered.qza \
  --o-filter-stats demux-joined-filter-stats.qza
qiime metadata tabulate --m-input-file demux-joined-filter-stats.qza  \
  --o-visualization demux-joined-filter-stats.qzv

#生产feature table  自带去除嵌合体
#deblur在denoising时需要输入整齐一样长度的序列，所以需要trim成相同的长度。这里需要用到前面join summary里我们记下来的join起来的质量有保障的序列大概有多长的数值。D#eblur的开发者们建议设置一个质量分数开始迅速下降的长度（recommend setting this value to a length where the median quality score begins to drop too low）。于#是本例中--p-trim-length为240。
qiime deblur denoise-16S \
  --i-demultiplexed-seqs demux-joined-filtered.qza \
  --p-trim-length 240 \
  --p-sample-stats \
  --p-jobs-to-start 2 \
  --p-min-reads 1 \
  --o-representative-sequences rep-seqs.qza \
  --o-table feature-table.qza \
  --o-stats deblur-stats.qza

#结果可视化
qiime feature-table summarize \
  --i-table feature-table.qza \
  --o-visualization feature-table.qzv

qiime deblur visualize-stats \
  --i-deblur-stats deblur-stats.qza \
  --o-visualization deblur-stats.qzv

### 方法4. q2 vsearch otu
echo vsearch start
date
cd $workdir
mkdir vsearch
cd vsearch
qiime vsearch dereplicate-sequences \
  --i-sequences $workdir/deblur/demux-joined-filtered.qza \
  --o-dereplicated-table table.qza \
  --o-dereplicated-sequences rep-seqs.qza
echo vsearch denovo start
date
#denovo pick otu
qiime vsearch cluster-features-de-novo \
  --i-table table.qza \
  --i-sequences rep-seqs.qza \
  --p-perc-identity 0.99 \
  --o-clustered-table table-dn-99.qza \
  --o-clustered-sequences rep-seqs-dn-99.qza
echo vsearch close start
date
#close reference 
qiime vsearch cluster-features-closed-reference \
  --i-table table.qza \
  --i-sequences rep-seqs.qza \
  --i-reference-sequences 85_otus.qza \
  --p-perc-identity 0.85 \
  --o-clustered-table table-cr-85.qza \
  --o-clustered-sequences rep-seqs-cr-85.qza \
  --o-unmatched-sequences unmatched-cr-85.qza
echo vsearch open start
date
#open referenced
qiime vsearch cluster-features-open-reference \
  --i-table table.qza \
  --i-sequences rep-seqs.qza \
  --i-reference-sequences 85_otus.qza \
  --p-perc-identity 0.85 \
  --o-clustered-table table-or-85.qza \
  --o-clustered-sequences rep-seqs-or-85.qza \
  --o-new-reference-sequences new-ref-seqs-or-85.qza
echo vsearch open end
date
#denovo 去除嵌合体
qiime vsearch uchime-denovo \
  --i-table atacama-table.qza \
  --i-sequences atacama-rep-seqs.qza \
  --output-dir uchime-dn-out

### 特征表和代表序列统计
qiime feature-table summarize \
  --i-table table_biom.qza \
  --o-visualization table.qzv \
  --m-sample-metadata-file $fastmap
qiime feature-table tabulate-seqs \
  --i-data repset-seqs.qza \
  --o-visualization rep-seqs.qzv
# 下载qzv在线查看，dada2只有1千多个ASV