RNA-seq中的基因表達量計算和表達差異分析

差異分析的步驟：
1）比對；
2） read count計算；
3） read count的歸一化；
4）差異表達分析；

背景知識：
1）比對：
普通比對： BWA，SOAP
開大GAP比對：Tophat（Bowtie2）；
2） Read count(多重比對的問題）：
丟棄
平均分配
利用Unique region估計並重新分配
表達量計算的本質
目標基因表達量相對參照系表達量的數值。
參照的本質：
（ 1）假設樣本間參照的訊號值應該是相同的；
（ 2）將樣本間參照的觀測值校正到同一水平；
（ 3）從參照的數值，校正並推算出其他觀測量的值。 
例如：Qpcr:目標基因表達量（迴圈數）相對看家基因表達量（迴圈數）；RNA-seq:目標基因的表達量（測序reads數），相對樣本RNA總表達量（總測序量的reads數），這是最常用的標準。
歸一化的原因及處理原則：
1）基因長度
2）測序量
3）樣本特異性（例如，細胞mRNA總量，汙染等）前兩者使用普通的RPKM演算法就可以良好解決，關鍵是第三個問題，涉及到不同的演算法處理。 
RNA-Seq歸一化演算法的意義：
基因表達量歸一化：在高通量測序過程中，樣品間在資料總量、基因長度、基因數目、高表達基因分佈甚至同一個基因的不同轉錄本分佈上存在差別。因此不能直接比較表達量，必須將資料進行歸一化處理。 
RNA-seq差異表達分析的一般原則
1）不同樣品的基因總表達量相似
2）上調差異表達與下調差異表達整體數量相似（上下調差異平衡）
3）在兩組樣品中不受處理效應影響的基因， 表達量應該是相近的（差異不顯著）。
4）看家基因可作為表達量評價依據（ 待定） 

不同的演算法比較：
以什麼數值來衡量表達量：RPKM、FPKM、TPM
以什麼作為參照標準：TMM（edgeR軟體）、De seq矯正
RPKM：是Reads Per Kilobase per Million mapped reads的縮寫，代表每百萬reads中來自於某基因每千鹼基長度的reads數。

本質：
1）以reads數為計算單位；
2）對基因長度（基因間的比較）和總資料量（樣本間的比較）做矯正；
RPKM的弊端
1）由於可變剪下，同一基因有效轉錄區域長度未必相同（這個一般情況下可以不考慮，瞭解一下：Cufflinks軟體考慮了這個問題）優化策略：外顯子或轉錄本水平的表達量分析。
2） 使用reads數計算基因表達量有輕微誤差（這裡暫不展開，主要了解一下定義）優化策略：FPKM或 TPM
3） mRNA的總量未必相等。 
RPKM的優化：FPKm
 F = Fragment，即測序片段數量。這些片段都是從完整的cDNA打碎而來的；
本質：以文庫中的片段數量為計算單位在Paired-end測序中，一個fragment就是兩條PE reads構成的片段。由於是PE比對，理論上比SE比對更可靠。 
RPKM的優化：TPM
T = Transcripts
本質：以轉錄本的條數為計算單位。使用轉錄本的條數（或者說：轉錄本的測序深度），代替reads數，在一定條件下定量更準，尤其樣本間表達基因總數差異很大的時候（例如，對照樣本有1萬個基因表達，另外處理組僅有4000個基因表達）。 
mRNA總量未必相等
mRNA總量不等——細胞本身不同
例如：活躍組織vs休眠的組織；癌細胞vs正常細胞
mRNA總量不等——汙染
例如：核糖體汙染外源RNA汙染 
解決方法——不同演算法比較
其中歸一化演算法介紹：
1）Total Count（TC）：總reads數矯正
2）Upper Quartile（UQ）：上四分之一分位數（總reads）
矯正
3）Median（Med）；中位數（總reads數）矯正
4）Quantile (Q)：基因晶片軟體limma中的校正演算法；
5）RPKM：總reads數，但引入了基因長度
6）幾何平均數：Deseq軟體中的演算法；
7）TMM：edgeR軟體中的演算法；
8）RPKM
邏輯1：不同位置數值的穩定性不同

四分位數quartile:將資料按從小到大排列，並分成四等分，這樣得到3個分割點，第一個分割點叫做lowerquartile，第二個叫Media，第三個叫Upper quartile
很顯然，極大值具有極大不穩定性，而且可能會顯著影
響總體之和（假設，我們之中有個馬雲，我們的總收入
有什麼變化？）
 所以，Upper quartile和Median的數值，比總表達量之
和更加穩定，更適合作為參照。
邏輯2：表達量居中的基因的表達量值，其數值應該是相似的。
DESeq與edgeR，預設情況下都使用這一的邏輯校正。（DESeq and edgeR Bioconductor packages） 
Deseq：異常高表達的基因，會顯著影響細胞中的總mRNA的數量。類似的，如果樣本中受到不同程度的外源RNA，如病毒、真菌等的汙染，也會顯著影響樣本總mRNA數，導致RPMK值的誤差。對於這樣的問題，Deseq嘗試對資料進行矯正（矯正因子），使表達量處於中間位置的基因表達量應該是基本相同的（即使用表達量處於中間的基因表達量值作為參照，而減少高表達基因的作用）。 
Deseq： 校正因子=樣本表達中位數/所有樣本表達量中位數：回答了一個關鍵的問題：Deseq不同差異比較組間，計算得到的表達量值不同。因
為樣本在變化，“所有樣本表達量的中位數”也在變動。RPKM：總表達量為參照
Deseq：中位數為參照 

TMM（edgeR）：與Deseq類似，在去除高表達基因和差異最大的基因後，TMM也是要找到一個加權係數，使剩餘的基因在被矯正後差異倍數可能小。TMM的加權係數是基於兩兩樣本比較後推算獲得的（也就是兩組樣本的比較，將產生與這次比較相關的加權係數）。然後將所有基因除以這個加權係數，從而保證大部分表達量居中的基因表達量最相似。 
不同RNA-seq表達量歸一化演算法的區別
Deseq類的校正演算法：理論上更加穩定；但不同批次的比較會得到不同的表達量值，不利於進行多處理組/批次資料的統一分析（例如，趨勢分析、共表達分析）校正會掩蓋一些問題（例如：樣本汙染）
RPKM類的演算法： 容易受異常高表達基因、外源汙染等的干擾；但也更容易從結果的異常中，發現潛在問題；得到的表達量值是恆定的，多處理組/批次的資料可以合併分析。折中的方法：使用RPKM類的演算法，但需要人工檢查資料是否
異常。備註： Deseq軟體也可以關閉校正的功能。

實際經驗總結
總之：從多方面考慮，RPKM類演算法，如果合理使用，依然是最優的。具體問題具體分析：在遇到問題的時候，找到問題的來源，從而給出解決方案（沒有完美的流程，只有最佳解決方案）