① 轉錄組時間序列數據處理
所謂時序分析 (time series analysis) 在 data science 中是非常重要的一個方向。對大多數商業行為而言如果能夠通過已有不同時間數據來進行預測就有可能大大提高自己的勝率。通常時間序列數據會包括趨勢部分和不規則部分, 我們需要做的就是剔除不規則部分然後找到趨勢所在,再進行預測。在預測過程中通常可以採用移動平均法、局部加權回歸法、指數平滑法和自回歸整合移動平均等方法。
生物學的時間相關數據本身預測屬性和商業數據相比要弱很多。一種是單一條件的純時間序列,主要看不同基因的表達模式,根據相似的表達譜將基因歸為多個類有助於找到功能相似的基因。另一種情況是含有對照和處理的時間序列,需要再考察不同條件的差異基因。
關於時間序列轉錄組數據分析的工具,近三年來有兩篇偏綜述和測評類的文章(一個人寫的)。
在這兩篇文章中還是提到了一些工具,但其中有一些用到matlab(這軟體貴啊),有一些年久失修或者不維護或者和最新R版本不兼容,篩篩撿撿能用的且文章里認為還不錯的也就剩下三四個。
來自於 DESeq 的方法,下文中提到的 ImpluseDE2 和 MaSigPro 都使用了這種模型。
來自於 maSigPro 方法,所謂多項式回歸區別常見的線性回歸,會把一次特徵轉換成高次特徵的線性組合多項式,比使用直線擬合更加准確。但是到底用幾次方需要具體分析,次數過高會出現過擬合。在能夠解釋自變數和因變數關系的前提下,次數應該是越低越好,這也算是奧卡姆剃刀原則吧。
所謂自回歸是統計上一種處理時間序列的方法,用 至 來預測本期 的表現並假設它們為線性關系。簡單說就是用自己來預測自己,因為是從回歸分析中線性回歸發展而來只是用x預測x,所以叫自回歸。
同樣是來自於 DESeq 的方法,下文中提到的 ImpluseDE2 和 MaSigPro 也都使用了這種方法。
似然比檢驗 (likelihood ratio test,LRT) 用於比較兩個模型的擬合優度進而確定哪個模型與樣本數據擬合的更好。其中一個是具有一定數量項的完整模型,另一個是刪掉完整模型中一部分項的簡化模型。LRT 檢驗中,自由度等於在簡化模型中減少的模型參數數目,LR近似符合卡方分布。一個相對復雜的模型與一個簡單模型比較,如果可以顯著地適合一個特定數據集,那麼這個復雜模型的附加參數就能夠用在以後的數據分析中。
為了測試多個時間點的任何差異,可以使用包含時間因子的設計和時間因子在簡化公式中被刪除的另一個設計。對於包括對照和實驗組的時間序列,可以使用包含條件因子,時間因子和兩者相互作用的公式。在這種情況下,使用具有不包含相互作用項的簡化模型的似然比檢驗將測試該條件是否在參考水平時間點(time 0)之後的任何時間點可以誘導基因表達的變化。
EBseq-HMM 採用的方法,來自於 BEseq。
這個軟體最早發表在2007年,相對老一些好在目前仍然在維護,其主要目的是給時序數據進行基於 模糊聚類演算法 的聚類。我們常見的聚類演算法可以分為嚴格聚類(hard clustering)和模糊聚類(Fuzzy clustering )(也叫做寬松聚類 soft clustering)。嚴格聚類會將一個基因只聚到一類中,kmeans 就屬於嚴格聚類。而模糊聚類允許同一數據屬於多個不同的類,其聚類結果是一個數據對聚類中心的隸屬度,0到1之間。對於分類很開的數據使用嚴格聚類是沒問題的。但對於時序表達量數據來說,不同的類常常會有重疊,所以可以嘗試寬松聚類方法。演算法需要首先設定一些參數,若初始化參數不合適,可能影響聚類結果的正確性。
在使用 Mfuzz 時首先應該進行數據標准化處理 ,可以使用類似於 FPKM 或者 TPM 的表達結果也可以使用 DESeq2 矯正後的結果進行比較分析,另外不支持值為0的數據,所以需要加上 pseudocount 。除此之外,Mfuzz 接受的數據格式為 ExpressionSet,需要對矩陣進行轉換。
這個包只能進行聚類,是找不了有處理對照組的差異基因的。需要注意。
運行masigpro 主要有四步:
有兩點內容需要注意:對於無對照的單一時序數據處理方法;以及處理轉錄數據時的特殊參數。因為這個包不會對數據進行標准化,所以應該提前做好,使用 DESeq2 即可。
另外,在實際分析的時候可能會出現 glm.fit: algorithm did not converge 的警告。這是由於進行 logistic 回歸時,依照極大似然估計原則進行迭代求解回歸系數,glm 函數默認最大迭代次數是 25,當數據不太好時 25 次迭代可能還不收斂,一方面可以增大迭代次數。但當增大迭代次數仍然不收斂就需要對數據進行異常值檢驗等進一步處理。通常把一些表達量極低或者極高的基因刪除掉,這個問題就可以解決。
ImpulseDE2 是最近才出來的一個R包,在前面提到的綜述評測文章中認為這個包找時序數據中的差異基因效果最好,它可以用來解決兩類問題。
這個包中,有一個 plotHeatmap 函數,可以藉助 ComplexHeatmap 對數據整體進行熱圖的繪制同時提取不同類的基因,也可以使用 plotGenes 看某一個基因的表達情況。
在展示的熱圖中會出現四部分,包括 transient and transition trajectorie,其中每一種 tarjectorie 又包括 up 和 down 兩類。所謂的 transient 可以理解為時序數據在中間某一個時間點存在up 或者 down peak,即在某一個時間點存在表達的最大或者最小值;而所謂的 transient 可以理解為一個持續的變化,比如持續的升高或者持續的降低。
EBSeq-HMM 是基於 EBSeq 二次開發的工具,主要用於分析時序數據。在計算的時候首先基於負二項分布對參數進行估計,然後利用自回歸隱馬模型將基因的表達進行分類。比較神奇的是,最終給到的結果會標示為 Up-Up-Down-Down-Down 之類的若干 path,然後你可以選出你感興趣的 path 進行後續分析。
因為目前做的數據是沒有對照的單一時間序列數據,所以還不能體會哪一個找出的差異基因更准確些。但是如果只是想把所有的基因根據不同的時間點分為若干表達 pattern,似乎結合 Mfuzz 和 ImpulseDE2 就可以了。
當然,涉及到聚類,尤其是非監督聚類的時候通常主觀因素還是較強,如果能對關鍵基因或者數據有一個大致的估計預判操作起來會相對輕鬆些,如果沒有,可能就需要結合不同類的生物學意義等角度來找合適的聚類數目了。
http://a-little-book-of-r-for-time-series.readthedocs.io/en/latest/src/timeseries.html
http://a-little-book-of-r-for-biomedical-statistics.readthedocs.io/en/latest/
https://laranikalranalytics.blogspot.com/2018/07/time-series-analysis-with-documentation.html
https://www.displayr.com/smoothing-time-series-data/?utm_medium=Feed&utm_source=Syndication
https://www.analyticsvidhya.com/blog/2015/12/complete-tutorial-time-series-modeling/
② 如何進行轉錄組數據分析 小小知識站
首先您做的事晶元呢?還是轉錄組測序呢? 晶元我不是太了解,可能是封閉系統,目前看對於發現新轉錄本不是很有利。 若是做轉錄組測序,首先去除 序列,然後片段重疊。然後將將每一個read 到基因組,取得GO值,功能注釋,轉錄水平評估,功能富集及pathway分析,若對新發現的轉錄本感興趣還可以做轉錄本的功能預測及細胞 。希望有高手來評價---說的不對的盡管指出。
③ 轉錄組測序所得的數據,怎樣傳到NCBI的SRA資料庫
xshell連接NCBI伺服器 輸入命令 ftp
連接成功後輸入賬號 密碼
mkdir命令創建一個文件夾,以你的BioProject accession命名
FTP 連接成功後連接成功後
進入剛才建立的BioProject文件夾
把原始數據拖進去即可
傳輸過程中,如果一定時間(五分鍾?)沒有激活傳輸界面,可能會導致連接斷開
進度可能會一直停在99%,重新連接成功後,也一直卡在99%,不過似乎數據還是上傳成功了的
傳輸成功的會顯示linked,同時FTP里對應的文件會消失 傳輸成功的會顯示linked,同時FTP里對應的文件會消失
全部文件夾都傳輸成功了會變成queued
然後就等NCBI處理了