
本文從QIIME2的官網(wǎng)進(jìn)行查看翻譯
https://docs.qiime2.org/2021.8/tutorials/moving-pictures/
首先使用Keemei對(duì)于metadata進(jìn)行check。
輸入的數(shù)據(jù)格式應(yīng)為*.qza
0 工作流程

總體概念圖:

1 對(duì)得到的數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析
例子數(shù)據(jù)
wget https://docs.qiime2.org/2021.8/data/tutorials/moving-pictures/demux.qza
以下代碼可以讓你知道每個(gè)樣本有多少sequences, 也會(huì)統(tǒng)計(jì)每個(gè)位置的質(zhì)量分布
qiime demux summarize --i-data demux.qza --o-visualization demux.qzv
所有的QIIME2 可視化文件都會(huì)是*.qzv格式,需要使用qiime tools view來(lái)查看
qiime tools view demux.qzv
結(jié)果:


2 sequence質(zhì)量控制和特征表構(gòu)建
DADA2適用于Illumina測(cè)序平臺(tái)數(shù)據(jù)的測(cè)試和校正pipeline。
此質(zhì)量控制過(guò)程將額外過(guò)濾在測(cè)序數(shù)據(jù)中識(shí)別的任何 phiX 讀數(shù)(通常存在于標(biāo)記基因 Illumina 序列數(shù)據(jù)中),并將過(guò)濾嵌合序列。
2.1 sequence質(zhì)量控制
data2 denosie-single法需要使用兩個(gè)參數(shù):
1: –p-trim-left 截取左端低質(zhì)量序列,我們看上圖中箱線圖,左端質(zhì)量都很高,無(wú)低質(zhì)量區(qū),設(shè)置為0;
2 : –p-trunc-len 序列截取長(zhǎng)度,也是為了去除右端低質(zhì)量序列,我們需要去除質(zhì)量下不好的端,可以看到120以后,甚至中位數(shù)都下降至20以下,需要全部去
請(qǐng)將下圖的所有輸出文件名字去掉“-dada2”,否則后續(xù)代碼會(huì)出現(xiàn)很多錯(cuò)誤。

m與n替換為實(shí)際選擇值
如果忘記去除“-dada2”,再運(yùn)行下述命令,使用
mv rep-seqs-dada2.qza rep-seqs.qza
mv table-dada2.qza table.qza
輸出文件:

對(duì)DADA2質(zhì)控的結(jié)果進(jìn)行可視化的文件輸出:

將輸出:

2.2 特征構(gòu)建和匯總
質(zhì)量過(guò)濾步驟完成后,您將需要瀏覽生成的數(shù)據(jù)。并且決定下一步Sampling Depth 的值,即決定所有的樣本的長(zhǎng)度,如果低于這個(gè)長(zhǎng)度的將會(huì)被去除
tabulate-seqs
使用命令:
feature-table summarize #每個(gè)樣本有多少sequences,每個(gè)的分布和一個(gè)統(tǒng)計(jì);該命令將為您提供與每個(gè)樣本和每個(gè)特征相關(guān)聯(lián)的序列數(shù)、這些分布的直方圖以及一些相關(guān)的匯總統(tǒng)計(jì)信息。
feature-table tabulate-seqs #給出每個(gè)IDs與每個(gè)簇代表sequences映射,這樣更容易進(jìn)行blast。并提供鏈接以輕松地針對(duì) NCBI nt 數(shù)據(jù)庫(kù)對(duì)每個(gè)序列進(jìn)行 BLAST
應(yīng)該把下面的table.qza改為上述的table-data2.qza, sample-metadata.tsv是來(lái)自上述stats-dada2.qzv文件的保存數(shù)據(jù)(保存后,可以改為這個(gè)名字)

輸出文件:

3 生成用于系統(tǒng)發(fā)育多樣性分析的樹
QIIME 支持多種系統(tǒng)發(fā)育多樣性指標(biāo),包括 Faith 的系統(tǒng)發(fā)育多樣性以及加權(quán)和未加權(quán)的 UniFrac。除了每個(gè)樣本的特征計(jì)數(shù)(即FeatureTable[Frequency]QIIME 2 工件中的數(shù)據(jù))之外,這些指標(biāo)還需要一個(gè)將特征相互關(guān)聯(lián)的有根系統(tǒng)發(fā)育樹。此信息將存儲(chǔ)在Phylogeny[Rooted]QIIME 2 工件中。為了生成系統(tǒng)發(fā)育樹,我們將使用插件中的align-to-tree-mafft-fasttree管道q2-phylogeny。
首先,管道使用mafft程序?qū)ξ覀冎械男蛄袌?zhí)行多序列比對(duì),F(xiàn)eatureData[Sequence]以創(chuàng)建FeatureData[AlignedSequence]QIIME 2 工件。接下來(lái),管道屏蔽(或過(guò)濾)對(duì)齊以移除高度可變的位置。這些位置通常被認(rèn)為會(huì)給生成的系統(tǒng)發(fā)育樹增加噪音。之后,管道應(yīng)用 FastTree 從掩碼對(duì)齊生成系統(tǒng)發(fā)育樹。FastTree 程序創(chuàng)建了一個(gè)無(wú)根樹,因此在本節(jié)的最后一步中,中點(diǎn)生根用于將樹的根放置在無(wú)根樹中最長(zhǎng)的尖端到尖端距離的中點(diǎn)處。
qiime 支持多種多樣性分析的指標(biāo):Faith’s Phylogenetic Diversity and weighted and unweighted UniFrac。 FeatureTable[Frequency] 生成有根的樹,以align-to-tree-mafft-fasttree進(jìn)行
使用mafft 命令, 再使用FastTree 生產(chǎn)無(wú)根樹

生成:

4 Aplha和beta多樣性分析
多樣性是使用q2-diversity進(jìn)行
首先使用core-metrics-phylogenetic法計(jì)算alpha和beta的多樣性指標(biāo)

命令:d是需要選擇使用的樣本測(cè)序深度,推薦查看上述的table.qzv文件

5 查看每個(gè)樣本的微生物組成
5.1 Aplha稀疏ploting
我們將使用 qiime diversity alpha-rarefaction探索 alpha 多樣性作為采樣深度的函數(shù)。 此可視化工具在多個(gè)采樣深度計(jì)算一個(gè)或多個(gè) alpha 多樣性指標(biāo),步長(zhǎng)介于 1(可選地使用 --pmin-depth 控制)和作為 --p-max-depth 提供的值之間。 在每個(gè)采樣深度步驟,將生成 10 個(gè)稀疏表,并且將為表中的所有樣本計(jì)算多樣性度量。 可以使用 --piterations 控制迭代次數(shù)(在每個(gè)采樣深度計(jì)算的精簡(jiǎn)表)。 將在每個(gè)偶數(shù)采樣深度為每個(gè)樣本繪制平均多樣性值,并且如果meta與 --m-metadata-file 參數(shù)一起提供,則樣本可以根據(jù)結(jié)果可視化中的元數(shù)據(jù)進(jìn)行分組。

(有問(wèn)題,沒(méi)做出來(lái)),直接下載的官網(wǎng)結(jié)果文件,可視化:


可視化將有兩個(gè)圖。頂部圖是 alpha 稀疏圖,主要用于確定樣本的豐富度是否已被完全觀察或排序。如果圖中的線在沿 x 軸的某個(gè)采樣深度處看起來(lái)“變平”(即接近零的斜率),則表明收集超出該采樣深度的其他序列不太可能導(dǎo)致觀察的附加功能。如果圖中的線條沒(méi)有變平,這可能是因?yàn)檫€沒(méi)有完全觀察到樣本的豐富度(因?yàn)槭占男蛄刑伲蛘呖赡鼙砻魅匀淮嬖诖罅繙y(cè)序錯(cuò)誤在數(shù)據(jù)中(這被誤認(rèn)為是新穎的多樣性)。
當(dāng)按元數(shù)據(jù)對(duì)樣本進(jìn)行分組時(shí),此可視化中的底部圖很重要。它說(shuō)明了當(dāng)特征表被稀疏到每個(gè)采樣深度時(shí),每個(gè)組中剩余的樣本數(shù)。如果給定的采樣深度d大于樣本的總頻率s(即,為 sample 獲得的序列數(shù)s),則不可能s在采樣深度為 sample 計(jì)算多樣性度量d。如果一組中的許多樣本的總頻率低于d,則該組的平均多樣性在d頂部圖將不可靠,因?yàn)樗窃谙鄬?duì)較少的樣本上計(jì)算的。因此,在按元數(shù)據(jù)對(duì)樣本進(jìn)行分組時(shí),必須查看底部圖以確保頂部圖中顯示的數(shù)據(jù)可靠。
注意:
提供的值--p-max-depth應(yīng)通過(guò)查看table.qzv上面創(chuàng)建的文件中顯示的“每個(gè)樣本的頻率”信息來(lái)確定。一般而言,選擇一個(gè)位于中值頻率附近的值似乎效果很好,但如果生成的稀疏圖中的線似乎沒(méi)有變平,您可能希望增加該值,或者如果您看起來(lái)像減少該值由于低總頻率比最大采樣深度更接近最小采樣深度而丟失許多樣本。
6 分類分析(Taxonomoc analysis)
在接下來(lái)的部分中,我們將開始探索樣本的分類組成,并再次將其與樣本元數(shù)據(jù)聯(lián)系起來(lái)。此過(guò)程的第一步是為FeatureData[Sequence]QIIME 2 工件中的序列分配分類法。我們將使用預(yù)訓(xùn)練的樸素貝葉斯分類器和q2-feature-classifier插件來(lái)做到這一點(diǎn)。該分類器在 Greengenes 13_8 99% OTU 上進(jìn)行了訓(xùn)練,其中序列已被修剪為僅包含來(lái)自本分析中測(cè)序的 16S 區(qū)域(V4 區(qū)域,由 515F/806R 引物對(duì)結(jié)合)的 250 個(gè)堿基。我們將這個(gè)分類器應(yīng)用到我們的序列中,我們可以生成從序列到分類法的映射結(jié)果的可視化。
注意:
分類分類器在根據(jù)您的特定樣品制備和測(cè)序參數(shù)(包括用于擴(kuò)增的引物和序列讀數(shù)的長(zhǎng)度)進(jìn)行訓(xùn)練時(shí)表現(xiàn)最佳。因此,通常您應(yīng)該按照使用 q2-feature-classifier 訓(xùn)練特征分類器中的說(shuō)明來(lái)訓(xùn)練您自己的分類分類器。我們?cè)?a target="_blank">數(shù)據(jù)資源頁(yè)面上提供了一些常用分類器,包括基于 Silva 的 16S 分類器,但將來(lái)作者可能會(huì)停止提供這些分類器,以便讓用戶訓(xùn)練自己的分類器,這將與他們的序列數(shù)據(jù)最相關(guān)。