流程
image.png
1.安裝與啟動
source /share/disk5/lianm/basic_tool/anaconda3/bin/activate qiime2-2019.7
source ~/anaconda3/bin/activate ~/anaconda3/envs/qiime2-2019.4/
conda deactivate
2.文件準(zhǔn)備
（1）Metadata files
https://docs.qiime2.org/2019.7/tutorials/metadata/
metadata
（2）MANIFEST file
Manifest
（3）FASTQ file
3.拆分樣品
4.可視化數(shù)據(jù)質(zhì)量
qiime demux summarize --i-data /share/disk5/zhuqh/16S/QIIME2_190samples/data/single-end-demux_test.qza --o-visualization /share/disk5/zhuqh/16S/QIIME2_190samples/data/single-end-demux_test.qzv
5.DADA2（降噪與聚類）
背景：
DADA2濾除有噪聲的序列，校正不確定序列中的錯誤，去除嵌合序列，去除單體(singletons，出現(xiàn)頻率僅有一次的序列)，然后對這些序列進(jìn)行去冗余。
（1）嵌合體序列由來自兩條或者多條模板鏈的序列組成，在16S/18S/ITS 擴(kuò)增子測序的分析中，系統(tǒng)相似度極高，嵌合體可達(dá)1%-20%，需要去除嵌合體序列。示意圖如下：
image.png
參數(shù)：
--p-trunc-len-f，表示位置前面的序列將被截?cái)啵?br>
--p-trunclen-r，指示讀取的位置后面序列截?cái)啵?br>
--p-max-ee，之前序列中超過預(yù)期最大錯誤率將被丟棄（默認(rèn)值為2）；
--p-truncq，截?cái)嗟谝粋€(gè)位置質(zhì)量分?jǐn)?shù)等于或小于提供值的序列（默認(rèn)值為2）；
--ptrim-left-f and--p-trim-left-r，如果引物存在于輸入序列文件中，可選參數(shù)可以設(shè)置為引物序列的長度，以便去噪。
#Dada2去噪
nohup qiime dada2 denoise-single --i-demultiplexed-seqs /share/disk5/zhuqh/16S/QIIME2_190samples/data/single-end-demux_test.qza --p-trim-left 0 --p-trunc-len 0 --o-representative-sequences 0_rep-seqs.qza --o-table 0_table.qza --o-denoising-stats 0_denoising-single-end-demux_stat.qza --p-n-threads 20 --p-max-ee 8 &
#分析數(shù)據(jù)統(tǒng)計(jì)
qiime metadata tabulate --m-input-file 390_denoising-single-end-demux_stat.qza --o-visualization 390_denoising-single-end-demux_stat.qzv
#table表統(tǒng)計(jì)
qiime feature-table summarize --i-table table.qza --o-visualization table.qzv --m-sample-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv
#seq表統(tǒng)計(jì)
qiime feature-table tabulate-seqs --i-data rep-seqs.qza --o-visualization rep-seqs.qzv
結(jié)果：
去噪過程輸出兩個(gè)工件：一個(gè)表文件以及代表性的序列文件。
算法：
（1）首先將每個(gè)reads全部看作單獨(dú)的單元，Sequence相同的reads被納入一個(gè)sequence，reads個(gè)數(shù)即成為該sequence的豐度（abundance）（其實(shí)就是去冗余的過程）
（2）計(jì)算每個(gè)sequence豐度的p-value。當(dāng)最小的p-value低于設(shè)定的閾值時(shí)，將產(chǎn)生一個(gè)新partition。每一個(gè)sequence將會被歸入最可能生成該
sequence的partition。
（3）依次類推，完成分割歸并。
6.物種分類
feature-classifier包括三種不同的分類方法。
（1）classify-consensus-blast和classify-consensus-vsearch都是基于比對的方法，可以在N個(gè)最好的比對結(jié)果中找一致最高的用于分類。這些方法直接參考數(shù)據(jù)庫FeatureData[Taxonomy]和FeatureData[Sequence]文件，不需要預(yù)先訓(xùn)練。
（2）基于機(jī)器學(xué)習(xí)的分類方法是通過classify-sklearn實(shí)現(xiàn)，理論上可以應(yīng)用任何分類方法。必須訓(xùn)練這些分類器，例如，為了解哪些特征可以最好地區(qū)分每個(gè)分類學(xué)組，在分類過程中添加額外的步驟。分類器訓(xùn)練過程是參考數(shù)據(jù)庫和特異的標(biāo)記基因，和每個(gè)標(biāo)記基因/參考數(shù)據(jù)庫組合計(jì)算一次；然后該分類器可以多次使用而不需要重新訓(xùn)練！
#下載物種注釋（一個(gè)預(yù)先訓(xùn)練好的物種注釋分類器）
wget -O "gg-13-8-99-515-806-nb-classifier.qza" "https://data.qiime2.org/2019.4/common/gg-13-8-99-515-806-nb-classifier.qza"
# 物種分類（使用機(jī)器學(xué)習(xí)分類器為序列分配可能的物種注釋）
nohup qiime feature-classifier classify-sklearn --i-classifier gg-13-8-99-515-806-nb-classifier.qza --i-reads 0_rep-seqs.qza --o-classification 0_taxa/0_taxonomy.qza &
# 物種結(jié)果轉(zhuǎn)換表格，可用于查看feature-物種
qiime metadata tabulate --m-input-file /share/disk5/zhengx/suke_qiime/0_taxa/0_taxonomy.qza --o-visualization /share/disk5/zhengx/suke_qiime/0_taxa/0_taxonomy.qzv
#物種分類柱狀圖
qiime taxa barplot --i-table 0_table.qza --i-taxonomy /share/disk5/zhengx/suke_qiime/0_taxa/0_taxonomy.qza --m-metadata-file mapping.tsv --o-visualization /share/disk5/zhengx/suke_qiime/0_taxa/0_taxa-bar-plots.qzv
7.構(gòu)建系統(tǒng)進(jìn)化樹
進(jìn)化樹是基于多序列比對的結(jié)果，可展示豐富的信息，用于Alpha, Beta多樣性分析的輸入文件。
#多序列比對，將去噪序列與mafft對齊
qiime alignment mafft --i-sequences rep-seqs.qza --o-alignment aligned-rep-seqs.qza
#移除高變區(qū)
qiime alignment mask --i-alignment aligned-rep-seqs.qza --o-masked-alignment masked-aligned-rep-seqs.qza
#使用the FastTree method建樹
nohup qiime phylogeny fasttree --i-alignment 0_masked-aligned-rep-seqs.qza --o-tree 0_unrooted-tree.qza --p-n-threads 10 &
#無根樹轉(zhuǎn)換為有根樹（在中間點(diǎn)對樹進(jìn)行根化）
qiime phylogeny midpoint-root --i-tree unrooted-tree.qza --o-rooted-tree rooted-tree.qza
8.Alpha多樣性
Alpha多樣性是計(jì)算樣品內(nèi)物種組成，包括數(shù)量和豐度兩維信息。Alpha多樣性計(jì)算前需要對OTU表進(jìn)行標(biāo)準(zhǔn)化，因?yàn)椴煌瑴y序深度，檢測到的物種數(shù)量會不同。我們將OTU表重抽樣至相同數(shù)據(jù)量，以公平比較各樣品的物種數(shù)量。
#計(jì)算多樣性（生成一系列的系統(tǒng)發(fā)育和非系統(tǒng)發(fā)育多樣性度量；為了比較序列深度不均勻的樣本，樣本中的最小序列數(shù)可用作次采樣深度，也可以稍微降低一點(diǎn)，即根據(jù)feature表的統(tǒng)計(jì)結(jié)果中Frequency per sample中minimum frequency，另外輸出文件夾需要提取mkdir core-metrics-results）
nohup qiime diversity core-metrics-phylogenetic --i-phylogeny 0_rooted-tree.qza --i-table 0_table.qza --p-sampling-depth 22663 --m-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv --output-dir 0_core-metrics-results &
# 統(tǒng)計(jì)Alpha多樣性的Faith’s phylogenetic diversity組間差異是否顯著，輸入多樣性值、實(shí)驗(yàn)設(shè)計(jì)，輸出統(tǒng)計(jì)結(jié)果。原理是所有組別和成對Kruskal Wallis檢驗(yàn)，一個(gè)非參數(shù)方差分析。
nohup qiime diversity alpha-group-significance --i-alpha-diversity 0_core-metrics-results/faith_pd_vector.qza --m-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv --o-visualization 0_core-metrics-results/faith-pd-group-significance.qzv &
# 統(tǒng)計(jì)evenness組間差異是否顯著
nohup qiime diversity alpha-group-significance --i-alpha-diversity 0_core-metrics-results/evenness_vector.qza --m-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv --o-visualization 0_core-metrics-results/sample-metadata.tsv &
#稀疏曲線：反應(yīng)數(shù)據(jù)的飽和度以及組內(nèi)多樣性，一看根據(jù)曲線逐漸持平測試測序深度是否足夠，二是最上面曲線的類別系統(tǒng)發(fā)育多樣性明顯高于其他類別；
qiime diversity alpha-rarefaction --i-table 0_table.qza --i-phylogeny 0_rooted-tree.qza --p-max-depth 80000 --m-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv --o-visualization 0_alpha-rarefaction.qzv
9.Beta多樣性
Beta多樣性是計(jì)算各樣品間的相同或不同.
#不同bodysite的unweighted unifrac距離的beta diversity 差異分析
qiime diversity beta-group-significance --i-distance-matrix 0_2_core-metrics-results/unweighted_unifrac_distance_matrix.qza --m-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv --m-metadata-column Group --o-visualization 0_2_core-metrics-results/unweighted-unifrac-body-site-significance.qzv --p-pairwise
#統(tǒng)計(jì)beta多樣性的組間差異是否顯著，PERMANOVA analysis on the BrayCurtis差異分析（通過PERMANOVA和ANOSIM統(tǒng)計(jì)方法）
qiime diversity beta-group-significance --i-distance-matrix core-metrics-results/unweighted_unifrac_distance_matrix.qza --m-metadata-file sample-metadata.tsv --m-metadata-column Subject --o-visualization core-metrics-results/unweighted-unifrac-subject-group-significance.qzv --p-pairwise
#基于unweighted-unifrac距離的beta diversity圖形可視化
qiime emperor plot --i-pcoa 0_2_core-metrics-results/unweighted_unifrac_pcoa_results.qza --m-metadata-file /share/disk5/zhuqh/16S/3groups_with_clinical_info_190/mapping.tsv --p-custom-axes DaysSinceExperimentStart --o-visualization 0_2_core-metrics-results/unweighted-unifrac-emperor-DaysSinceExperimentStart.qzv
#基于bray-curtis距離的beta diversity圖形可視化
qiime emperor plot --i-pcoa core-metrics-results/bray_curtis_pcoa_results.qza --m-metadata-file sample-metadata.tsv --p-custom-axes DaysSinceExperimentStart --o-visualization core-metrics-results/bray-curtis-emperor-DaysSinceExperimentStart.qzv
【附錄思路】
看測序儀，illumine是正常掉峰，其余測序儀質(zhì)量結(jié)果不規(guī)范；以及保存的資料；
去掉數(shù)據(jù)只影響豐度；
？高峰度的OTU？基于注釋結(jié)果聚類，到門或
？a多樣性看飽和度，是否分析全
？如何挖信息（個(gè)體間差異大，臨床信息，注釋文件是啥，根據(jù)信息重新分組）
？臨床基本信息，數(shù)據(jù)量，導(dǎo)致尿感是一種還是多種？原因是什么？尿培養(yǎng)和16s結(jié)果對比（一致和不一致的比例？注釋的情況？根據(jù)物種可以把原參與人群分為幾類；
mkdir temp
source /software/miniconda/activate /software/.../qiime2-2019.4/
qiime tools export --input-path table.qza --output-path ./temp/
cd out
#格式轉(zhuǎn)換
biom convert -i feature-table.biom -o oyold2.txt
2.qimme2上課所有命令
#請大家在D盤下面新建文件夾qiime2，把data文件夾、gg-13-8-99-515-806-nb-classifier.qza以及sample-metadata.tsv放在qiime2文件夾中
#輸入文件準(zhǔn)備
cd share
ll data
#加載數(shù)據(jù)
qiime tools import --type EMPSingleEndSequences --input-path data --output-path emp-single-end-sequences.qza
#Demultiplexing sequences
qiime demux emp-single --i-seqs emp-single-end-sequences.qza --m-barcodes-file sample-metadata.tsv --m-barcodes-column BarcodeSequence --o-per-sample-sequences demux.qza
qiime demux summarize --i-data demux.qza --o-visualization demux.qzv
qiime tools view demux.qzv
#用DADA2方法進(jìn)行質(zhì)控（耗時(shí)）
qiime dada2 denoise-single --i-demultiplexed-seqs demux.qza --p-trim-left 0 --p-trunc-len 120 --o-representative-sequences rep-seqs.qza --o-table table.qza
#生成feature文件
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
qiime tools view table.qzv
#phylogenetic diversity
#多序列比對
qiime alignment mafft --i-sequences rep-seqs.qza --o-alignment aligned-rep-seqs.qza
#移除多變區(qū)
qiime alignment mask --i-alignment aligned-rep-seqs.qza --o-masked-alignment masked-aligned-rep-seqs.qza
#構(gòu)建無根樹
qiime phylogeny fasttree --i-alignment masked-aligned-rep-seqs.qza --o-tree unrooted-tree.qza
#無根樹轉(zhuǎn)為有根樹
qiime phylogeny midpoint-root --i-tree unrooted-tree.qza --o-rooted-tree rooted-tree.qza
#Alpha and beta diversity analysis
#生成alpha diversity不同指數(shù)和beta diversity不同距離下的結(jié)果：sampling depth 一般為數(shù)據(jù)量最小的樣本的序列數(shù)
qiime diversity core-metrics-phylogenetic --i-phylogeny rooted-tree.qza --i-table table.qza --p-sampling-depth 1109 --m-metadata-file sample-metadata.tsv --output-dir core-metrics-results
#進(jìn)化多樣性差異分析
qiime diversity alpha-group-significance --i-alpha-diversity core-metrics-results/faith_pd_vector.qza --m-metadata-file sample-metadata.tsv --o-visualization core-metrics-results/faith-pd-group-significance.qzv
qiime tools view core-metrics-results/faith-pd-group-significance.qz
#均勻度差異分析
qiime diversity alpha-group-significance --i-alpha-diversity core-metrics-results/evenness_vector.qza --m-metadata-file sample-metadata.tsv --o-visualization core-metrics-results/evenness-group-significance.qzv
#不同bodysite的unweighted unifrac距離的beta diversity 差異分析
qiime diversity beta-group-significance --i-distance-matrix core-metrics-results/unweighted_unifrac_distance_matrix.qza --m-metadata-file sample-metadata.tsv --m-metadata-column BodySite --o-visualization core-metrics-results/unweighted-unifrac-body-site-significance.qzv --p-pairwise
#不同subject的unweighted unifrac距離的beta diversity 差異分析
qiime diversity beta-group-significance --i-distance-matrix core-metrics-results/unweighted_unifrac_distance_matrix.qza --m-metadata-file sample-metadata.tsv --m-metadata-column Subject --o-visualization core-metrics-results/unweighted-unifrac-subject-group-significance.qzv --p-pairwise
#基于unweighted-unifrac距離的beta diversity圖形可視化
qiime emperor plot --i-pcoa core-metrics-results/unweighted_unifrac_pcoa_results.qza --m-metadata-file sample-metadata.tsv --p-custom-axes DaysSinceExperimentStart --o-visualization core-metrics-results/unweighted-unifrac-emperor-DaysSinceExperimentStart.qzv
#基于bray-curtis距離的beta diversity圖形可視化
qiime emperor plot --i-pcoa core-metrics-results/bray_curtis_pcoa_results.qza --m-metadata-file sample-metadata.tsv --p-custom-axes DaysSinceExperimentStart --o-visualization core-metrics-results/bray-curtis-emperor-DaysSinceExperimentStart.qzv
#Alpha rarefaction plotting
#稀疏曲線：反應(yīng)數(shù)據(jù)的飽和度以及組內(nèi)多樣性
qiime diversity alpha-rarefaction --i-table table.qza --i-phylogeny rooted-tree.qza --p-max-depth 4000 --m-metadata-file sample-metadata.tsv --o-visualization alpha-rarefaction.qzv
#Taxonomic analysis（耗時(shí)）
qiime feature-classifier classify-sklearn --i-classifier gg-13-8-99-515-806-nb-classifier.qza --i-reads rep-seqs.qza --o-classification taxonomy.qza
qiime metadata tabulate --m-input-file taxonomy.qza --o-visualization taxonomy.qzv
qiime taxa barplot --i-table table.qza --i-taxonomy taxonomy.qza --m-metadata-file sample-metadata.tsv --o-visualization taxa-bar-plots.qzv
#ANCOM進(jìn)行差異分析
#對于gut這一組的樣本進(jìn)行差異分析
qiime feature-table filter-samples --i-table table.qza --m-metadata-file sample-metadata.tsv --p-where "BodySite='gut'" --o-filtered-table gut-table.qza
#去掉0值
qiime composition add-pseudocount --i-table gut-table.qza --o-composition-table comp-gut-table.qza
#不同subject差異分析
qiime composition ancom --i-table comp-gut-table.qza --m-metadata-file sample-metadata.tsv --m-metadata-column Subject --o-visualization ancom-Subject.qzv
#合并屬水平 對gut樣本進(jìn)行差異分析
qiime taxa collapse --i-table gut-table.qza --i-taxonomy taxonomy.qza --p-level 6 --o-collapsed-table gut-table-l6.qza
#去掉0值
qiime composition add-pseudocount --i-table gut-table-l6.qza --o-composition-table comp-gut-table-l6.qza
#不同subject中
qiime composition ancom --i-table comp-gut-table-l6.qza --m-metadata-file sample-metadata.tsv --m-metadata-column Subject --o-visualization l6-ancom-Subject.qzv
#文件輸出
#輸出feature table
qiime tools export table.qza --output-dir out
cd out
#biom格式轉(zhuǎn)換
biom convert -i feature-table.biom -o oyold2.txt --to-tsv
qiime metadata tabulate --m-input-file rep-seqs.qza --m-input-file taxonomy.qza --o-visualization tabulated-feature-metadata.qzv
qiime tools view tabulated-feature-metadata.qzv
