設(shè)置生物學重復這個環(huán)節(jié)也是你實驗設(shè)計很重要的一part，設(shè)置的好對你下游分析也有利，通常我們做轉(zhuǎn)錄組測序，需要的樣本量每組至少為3個生物學重復，這個處理起來就很合理，并且現(xiàn)在流行的差異分析軟件DEseq2，limma，edgeR等等都是針對有重復的數(shù)據(jù)去做的，但有時候會不幸碰到樣品測序失敗不能用，導致每組就給你剩一個重復時候該怎么辦，之前我有批數(shù)據(jù)就是這樣，但是辦法總比困難多不能放過任何實驗數(shù)據(jù)，搜了搜其實還是有一些方法可以去解決的，在這里介紹下我搜到的幾種方法。
假如現(xiàn)在你手頭有如下文件(test.txt)，只有倆樣品RPKM_A (對照) 和RPKM_B (處理), 值為標準化后的RPKM。
1. 根據(jù)foldchange直接篩選
之前在一篇中文文獻中見到有人用這種方法，作者自定義差異基因的標準：至少有一組RPKM值大于5，且滿足foldchange(差異倍數(shù)) > 2，我們可以在LInux中直接可以用awk進行過濾，其實Excel、R中也可以操作，根據(jù)個人習慣吧。代碼如下：
### 上調(diào)基因########
# 提取B組大于等于5，A組等于0的基因。
less test.txt | gawk '{if (($2==0)&&($3>=5)) print $0}' > up.txt
# 提取A、B倆組至少有一組大于等于5，且B組值/A組值大于等于2
less test.txt | gawk '{if (($2!=0)&&($3!=0)) print $0}'|gawk '{if (($2>=5)||($3>=5)) print $0}'|sed '1d'|gawk '{if ($3/$2>=2) print $0}' >> up.txt
### 下調(diào)基因#########
# 提取A組大于等于55，B組等于0的基因
less test.txt | gawk '{if (($2>=5)&&($3==0)) print $0}' > down.txt
# 提取A、B倆組至少有一組大于等于5，且A組值/B組值大于等于2
less test.txt | gawk '{if (($2!=0)&&($3!=0)) print $0}'|gawk '{if (($2>=5)||($3>=5)) print $0}'|sed '1d'|gawk '{if ($2/$3>=2) print $0}' >> down.txt
2. edgeR包
這種方法我在提到過，edgeR包可以做無重復的差異分析，不過需要認為指定一個dispersion值(設(shè)置BCV值)，這樣得到的結(jié)果比較主觀，不同的人就可以有不同的結(jié)果。通常如果是實驗控制的好的人類數(shù)據(jù)，那么選擇BCV=0.4，比較好的模式生物選擇BCV=0.1。參考
代碼如下：
library(edgeR)
##跟DESeq2一樣，導入數(shù)據(jù)，預處理（用了cpm函數(shù)）
exprSet<- read.table(file = "test.txt", sep = "\t", header = TRUE, row.names = 1, stringsAsFactors = FALSE)
group_list <- factor(c(rep("Contral",1),rep("Treat",1)))
##設(shè)置分組信息，并做TMM標準化
exprSet <- DGEList(counts = exprSet, group = group_list)
bcv = 0.1 #設(shè)置BCV值
et <- exactTest(exprSet, dispersion=bcv^2)
write.csv(topTags(et, n = nrow(exprSet$counts)), 'result.csv', quote = FALSE) #輸出主要結(jié)果
結(jié)果文件如下：
3. Gfold軟件
地址：https://zhanglab.#edu.cn/softwares/GFOLD/index.html
Gfold軟件應(yīng)該是做沒有生物學重復樣本首選的軟件，該軟件由同濟大學開發(fā)的，網(wǎng)站 往下拉可以看到該軟件的幾個功能，其中Example3為鑒定無重復的數(shù)據(jù)的差異基因。另外，這個軟件不支持Windows 版本，是基于linux的一個安裝軟件，所以我們需要在linux環(huán)境下使用
安裝GSL
使用Gfold之前必須先安裝Gsl,如圖下載最新的版本
# 安裝最新版的
wget http://mirrors.ocf.berkeley.edu/gnu/gsl/gsl-latest.tar.gz
tar -zxv -f gsl-latest.tar.gz
.configure --prefix=/home/pub_guest/hekai/soft_ware/GFOLD/gsl-2.6/
make
make check(選做)
make install
安裝Gfold
- 下載安裝包
wget https://zhanglab.#edu.cn/softwares/GFOLD/gfold.V1.1.4.tar.gz
tar -zxv - f gfold.V1.1.4.tar.gz
cd gfold.V1.1.4
- 我們打開目錄下的README文件，可以需要執(zhí)行下面?zhèn)z句。注意 /your/installed/path/ 是你安裝GSL的路徑，把下面命令行中的替換為相應(yīng)的路徑即可
# 打開你的bashrc文件，再最后添加
export CXXFLAGS="-g -O3 -I/your/installed/path/include -L/your/installed/path/lib"
export LD_LIBRARY_PATH="/your/installed/path/lib:"$LD_LIBRARY_PATH
source ~/.bashrc
- 接著輸入make, 發(fā)現(xiàn)報錯了。
- 不慌，文檔里也有提醒如果報錯 直接輸入以下命令行即可
# If it happens, follow step 1 again. If error remains, try the following command:
g++ -O3 -Wall -g main.cc -o gfold -lgsl -lgslcblas -I/your/installed/path/include -L/your/installed/path/lib
-
這是我們就會發(fā)現(xiàn)目錄下多了一個gfold軟件，是可執(zhí)行狀態(tài)，我們.gfold -h ，可以看到該軟件的幫助文檔，證明此時已經(jīng)安裝成功了，我們將其添加為環(huán)境變量里方便我們使用。
echo 'export PATH=/your/installed/path:$PATH' >>~/.bashrc ###### **/your/installed/path/** 是你安裝**gfold.V1.1.4**的路徑 source ~/.bashrc
差異分析
- 我們需要準備倆個文件Control.txt和Treat.txt，我們看下處理組Control.txt文件都包含哪些，Gfold輸入文件規(guī)定必須為5列，前倆列可以分別輸入你的基因ID號和Symbol號，倆列內(nèi)容一樣其實也不影響，第三列為原始Counts值，第四列為基因長度，最后一列為標準化的RPKM值，同樣對照組Treat.txt文件也按照這樣準備。
-
準備完畢后，直接一條命令計算差異。
gfold diff -s1 Control.txt -s2 Treat.txt -o ControlVSTreat.csvimage-20200928162344142 -
OK,已經(jīng)計算完了，我們看下結(jié)果文件都有哪些內(nèi)容。主要一共7列信息，前兩列沒什么可說，就是gene symbol和gene name，第三列是GFOLD值，相當于log2（Fold Change），該值等于0的基因則記為非差異基因，非0的值才是差異基因，E-FDR是基于重復的Empirical FDR，因此無重復樣本的經(jīng)驗FDR均為1。Log2fdc以及后面的RPKM列可以忽略考慮，因為最開始的exon的長度，我們是給定的是一個虛擬的數(shù)據(jù)。所以真正的確定差異是否顯著，主要是看GFOLD值，GFOLD>0,代表case組中高表達，GFOLD<0,代表case組中低表達。后續(xù)篩選差異基因如果太多或者太少，我們也可以通過設(shè)定GFOLD的閾值來控制，比如設(shè)定<-0.3或者大于0.3。
image-20200928162708344
除上述之外幾種方法，還有幾個比較經(jīng)典的軟件可供選擇進行無重復數(shù)據(jù)的差異分析，比如BMC Bioinformatics發(fā)表過的一篇文章中提出了一種新的差異基因分析方法利用貝葉斯方法來推斷真實基因表達數(shù)的 后驗分布。其創(chuàng)新型之一該方法包括了由RNA樣品濃度決定的覆蓋度參數(shù)，之二是真實基因表達量后驗分布的比較為尋找差異表達基因提供了一個參照。這種方法針對無重復樣本的數(shù)據(jù)是有一定優(yōu)勢的，這里提供一個鏈接大家有興趣的話可以去看該博主的講解,之后有機會也會嘗試一下該軟件的使用進行比較。
