梅洛斯
@swallow2038
分享
Mon, Feb 8, 2021 12:24 PM
Mon, Feb 8, 2021 2:51 PM
1
持續更新
[生資學習筆記]Day6
實作篇part.3
梅洛斯
@swallow2038
Mon, Feb 8, 2021 12:24 PM
前提提要
@swallow2038 - [生資學習筆記] Day4 實作篇PART.1 資料來自dada2官方網站
@swallow2038 - 持續更新 [生資學習筆記]Day5 實作篇PART.2 PART.1在這裡...
梅洛斯
@swallow2038
Mon, Feb 8, 2021 12:25 PM
今天預計要講我覺得最麻煩的章節: 物種比對
梅洛斯
@swallow2038
Mon, Feb 8, 2021 12:35 PM
先說說物種比對是做什麼用的
簡單講是拿已知資料庫的序列，跟自己定序組裝出來的序列做比對
之所以麻煩在於每種資料庫的數據可能長得都不一樣
即使是16S最常用的SILVA資料庫裏頭還有好幾種版本供使用者選用
梅洛斯
@swallow2038
Mon, Feb 8, 2021 12:39 PM
dada2選用的SILVA資料庫 是用機器學習的訓練集(training set)，裡面的資料經過驗證、測試視為一個好的模型才丟出來給大家使用
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:17 PM
前面數據經過篩選、合併、去除不應該有的chimera序列後，才會進行到物種比對這邊
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:19 PM
進入dada2網站提供的database下載區，點選最新v138 SILVA版本後再往下做
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:22 PM
參考資料庫像似標準答案，有既定已知物種序列資訊
如果比對成功的話會寫入資料至該序列
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:25 PM
指令版本1: assignTaxonomy + addSpecies
assignTaxonomy 顧名思義做物種比對，範例中的database只有到Genus屬的資料，Species種的資料要另外下載再執行addSpecies指令才會寫入對應資料。
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:28 PM
選用silva_nr99_v138_train_set.fa.gz 資料集進行運算
後續資料再以silva_species_assignment_v138.fa.gz
輸入Species資料
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:29 PM
指令版本2: assignTaxonomy
官方提供silva_nr99_v138_wSpecies_train_set.fa.gz
裡面含有Species資料集 (With Species)，所以理論上跑一次資料集就好
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:31 PM
除了這兩版本之外，還有其他套件如DECIPHER的資料集可以套用
指令為IdTaxa，也需要官方提供的數據集供比對
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:32 PM
每個套件裡面的training set我覺得都不一樣
即使都是SILVA資料庫做轉換的檔案，跑出來結果就是不一樣
像似可樂的配方一樣神秘
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:35 PM
跑完後只要檔案輸出成功，再把資料欄位轉換成ASVs
用write.table匯出你指定的名稱，這部分就大功告成了
後面才是麻煩的開始
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:36 PM
如果是既定答案的標準菌株
data還好解釋
但像土壤微生物或是腸內菌落種類百百種，不見得在database找得到，能鑑定到Genus屬階級已經是很詳細的資料了
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:42 PM
實質上知道Genus種就能代表樣品中真的存在某菌株嗎?
只能說推斷有，因為環境微生物8成以上是實驗室培養不出來的。
幾個major group做鑑定是沒有太大的問題，但即使是ASVs號稱能把PCR error降至最低，樣品中只有少數存在的物種可能被其他微生物"占領資源"，使得序列被判定有問題。
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:43 PM
還有不同database鑑定物種的名稱用法不一樣
有碰到分析完結果60%全都是NA (not assigned) = 無名
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:45 PM
因為沒有100%最好的pipeline
只有最適合自己家的pipeline
這個只能自己try troubleshooting才知道
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:50 PM
示意圖
這是我自己try三種版本的資料集跑完的結果
total ASVs大家都是固定一樣的
用俗稱Venn diagram的交集圖呈現時
左上跟右上兩者數值幾乎一致
但下方綠色圈圈的資料明顯有90幾筆跟別人沒有交集
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:53 PM
最開始我是用最下方的綠色IdTaxa跑數據
結果資料出現好幾筆"NA"對不起來
這跟演算法&訓練集有關係
而且IdTaxa的處理結果是沒有到Species種的資料
想偷吃步用別人家的fasta檔案還不行
像是安卓不能跟蘋果手機互換這樣
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:57 PM
總之
物種鑑定這邊運行&除錯花的時間比我想像中還久
而且問題基本上沒有制式的解答
如果你做5重複試驗，每個樣品結果有些物種有出現但有些不一樣
變成你要如何取捨data要保留還是刪除
這關係到後續數據重現能不能解釋你最開始試驗組處理的"效果"是否達到既定的目標
梅洛斯
@swallow2038
Mon, Feb 8, 2021 1:58 PM
嘛後半部變成心得&抱怨串了
先到這裡
載入新的回覆
[生資學習筆記]Day6
實作篇part.3
簡單講是拿已知資料庫的序列，跟自己定序組裝出來的序列做比對
之所以麻煩在於每種資料庫的數據可能長得都不一樣
即使是16S最常用的SILVA資料庫裏頭還有好幾種版本供使用者選用
進入dada2網站提供的database下載區，點選最新v138 SILVA版本後再往下做
如果比對成功的話會寫入資料至該序列
assignTaxonomy 顧名思義做物種比對，範例中的database只有到Genus屬的資料，Species種的資料要另外下載再執行addSpecies指令才會寫入對應資料。
後續資料再以silva_species_assignment_v138.fa.gz
輸入Species資料
官方提供silva_nr99_v138_wSpecies_train_set.fa.gz
裡面含有Species資料集 (With Species)，所以理論上跑一次資料集就好
指令為IdTaxa，也需要官方提供的數據集供比對
即使都是SILVA資料庫做轉換的檔案，跑出來結果就是不一樣
像似可樂的配方一樣神秘用write.table匯出你指定的名稱，這部分就大功告成了
後面才是麻煩的開始data還好解釋
但像土壤微生物或是腸內菌落種類百百種，不見得在database找得到，能鑑定到Genus屬階級已經是很詳細的資料了
只能說推斷有，因為環境微生物8成以上是實驗室培養不出來的。
幾個major group做鑑定是沒有太大的問題，但即使是ASVs號稱能把PCR error降至最低，樣品中只有少數存在的物種可能被其他微生物"占領資源"，使得序列被判定有問題。
有碰到分析完結果60%全都是NA (not assigned) = 無名
只有最適合自己家的pipeline
這個只能自己try troubleshooting才知道
這是我自己try三種版本的資料集跑完的結果
total ASVs大家都是固定一樣的
用俗稱Venn diagram的交集圖呈現時
左上跟右上兩者數值幾乎一致
但下方綠色圈圈的資料明顯有90幾筆跟別人沒有交集
結果資料出現好幾筆"NA"對不起來
這跟演算法&訓練集有關係
而且IdTaxa的處理結果是沒有到Species種的資料
想偷吃步用別人家的fasta檔案還不行
像是安卓不能跟蘋果手機互換這樣物種鑑定這邊運行&除錯花的時間比我想像中還久
而且問題基本上沒有制式的解答
如果你做5重複試驗，每個樣品結果有些物種有出現但有些不一樣
變成你要如何取捨data要保留還是刪除
這關係到後續數據重現能不能解釋你最開始試驗組處理的"效果"是否達到既定的目標
先到這裡