為什么Alpha多樣性的輸入數(shù)據(jù)會(huì)是它?
生物信息學(xué)習(xí)的正確姿勢(shì)
NGS系列文章包括NGS基礎(chǔ)、在線繪圖、轉(zhuǎn)錄組分析?(Nature重磅綜述|關(guān)于RNA-seq你想知道的全在這)、ChIP-seq分析?(ChIP-seq基本分析流程)、單細(xì)胞測(cè)序分析?(重磅綜述:三萬字長(zhǎng)文讀懂單細(xì)胞RNA測(cè)序分析的最佳實(shí)踐教程)、DNA甲基化分析、重測(cè)序分析、GEO數(shù)據(jù)挖掘(典型醫(yī)學(xué)設(shè)計(jì)實(shí)驗(yàn)GEO數(shù)據(jù)分析 (step-by-step))、批次效應(yīng)處理等內(nèi)容。
在擴(kuò)增子培訓(xùn)過程中,總會(huì)被問到一個(gè)問題:抽平后的OTU表和計(jì)算相對(duì)豐度后的OTU表差別是什么?
抽平和計(jì)算相對(duì)豐度都是對(duì)OTU表進(jìn)行標(biāo)準(zhǔn)化的方式,而抽平后的OTU表一般就會(huì)用于計(jì)算Alpha多樣性,兩者都可用于計(jì)算樣本間的差異OTU。但是為什么不統(tǒng)一使用抽平后的OTU表或計(jì)算相對(duì)豐度的OTU表呢?
我們?cè)?a target="_blank" textvalue="擴(kuò)增子培訓(xùn)" data-itemshowtype="11" tab="innerlink" data-linktype="2">擴(kuò)增子培訓(xùn)中學(xué)到了兩種以上的方式計(jì)算Alpha多樣性,比如用vegan包去計(jì)算6種Alpha多樣性、usearch計(jì)算14種Alpha多樣性等。下面我們可以嘗試使用抽平后的OTU表(otutab_rare.txt)和計(jì)算相對(duì)豐度后的OTU表(otutab.freq.txt)分別來計(jì)算Alpha多樣性。
方法一:vegan包計(jì)算Alpha多樣性
suppressWarnings(suppressMessages(library(vegan)))
otu_rare=read.table("otutab_rare.txt", header=T, sep="\t", quote = "", row.names=1, comment.char="")
otu_freq=read.table("otutab.freq.txt", header=T, sep="\t", quote = "", row.names=1, comment.char="")
# 計(jì)算抽平后的OTU表的richness, chao1和ACE指數(shù)
estimateR_rare = t(estimateR(t(otu_rare)))[,c(1,2,4)]
head(estimateR_rare)
## S.obs S.chao1 S.ACE
## KO1 2148 2442.130 2442.174
## KO2 2097 2376.372 2387.770
## KO3 1792 2169.230 2147.837
## KO4 1853 2174.955 2199.956
## KO5 1657 2061.449 2091.146
## KO6 1976 2303.736 2294.931
# 計(jì)算相對(duì)豐度后的OTU表的richness, chao1和ACE指數(shù)
## 由于代碼運(yùn)行報(bào)錯(cuò)所以就注釋掉了,報(bào)錯(cuò)內(nèi)容在注釋中進(jìn)行展示
#estimateR_freq = t(estimateR(t(otu_freq)))[,c(1,2,4)]
## Error in estimateR.default(newX[, i], ...) : function accepts only integers (counts)方法二:usearch計(jì)算Alpha多樣性
# input rare table to calculate the alpha diversty
usearch -alpha_div otutab_rare.txt -output rare.alpha.txt
## usearch v10.0.240_win32, 2.0Gb RAM (8.4Gb total), 8 cores
## (C) Copyright 2013-17 Robert C. Edgar, all rights reserved.
## http://drive5.com/usearch
##
## License: [email protected]
##
## 00:00 4.5Mb 0.1% Reading otutab_rare.txt00:00 5.2Mb 100.0% Reading otutab_rare.txt
# input otu rare table to calculate the alpha diversty
#usearch -alpha_div otutab.freq.txt -output freq.alpha.txt
## ---Fatal error---
##alphadiv.cpp(192) assert failed: m_TotalCount > 0可以看到當(dāng)輸入數(shù)據(jù)是抽平后的OTU表(otutab_rare.txt)時(shí),無論用哪種方法計(jì)算Alpha多樣性都沒有問題,而輸入數(shù)據(jù)是計(jì)算相對(duì)豐度后的OTU表(otutab.freq.txt)時(shí),代碼紛紛報(bào)錯(cuò)。
由于usearch并未開源,所以我們主要來看vegan的estimateR函數(shù)的報(bào)錯(cuò)原因function accepts only integers (counts),表明該函數(shù)只接受整數(shù),并且函數(shù)的文檔也指明Function estimateR is based on abundances (counts) on single sample site,說明輸入數(shù)據(jù)需要是豐度的counts值,但更具體的原因是什么呢?這就和Alpha多樣性的不同指標(biāo)的計(jì)算方式有關(guān)了。
estimateR函數(shù)計(jì)算了三種Alpha多樣性指標(biāo),分別是Richness、Chao1、ACE。
Richness是最好計(jì)算的Alpha多樣性指數(shù),其計(jì)算方式是

其中N就是該樣本的物種數(shù)目;
Chao1是常用的Alpha多樣性指數(shù)之一,其計(jì)算方式是

其中N是該樣本的物種數(shù)目;S是該樣本中豐度為1的物種數(shù)目;D是該樣本中豐度為2的物種數(shù)目;
ACE是常用的Alpha多樣性指數(shù)之一,其計(jì)算方式是


可以看到計(jì)算Chao1指數(shù)和ACE指數(shù)都依據(jù)豐度的counts值來計(jì)算,而輸入計(jì)算了相對(duì)豐度的OTU表時(shí)就丟失了豐度的counts值信息,故而只能使用抽平后的OTU表。
那你還知道哪些基于豐度的counts值計(jì)算的Alpha多樣性指標(biāo)或者Alpha多樣性其他容易被忽視又很重要的知識(shí)點(diǎn)呢?歡迎在評(píng)論區(qū)討論。

視頻演示Alpha多樣性圖形繪制
往期精品(點(diǎn)擊圖片直達(dá)文字對(duì)應(yīng)教程)
后臺(tái)回復(fù)“生信寶典福利第一波”或點(diǎn)擊閱讀原文獲取教程合集




























