簡單使用limma做差異分析
首先需要說明的是，limma是一個非常全面的用於分析芯片以及RNA-Seq的差異分析，按照其文章所說：
limma is an R/Bioconductor software package that provides an integrated solution for analysing data from gene expression experiments.
在這我只是對其中的一種情況進行簡單的總結，比如這個包可以處理RNA-Seq數據，我簡單的以兩個比較組進行分組為例，至於其他分組情況，請看limma說明文檔，有非常詳細的說明，非常親民。
-
首先我們還是輸入count矩陣，這里也跟其他差異分析R包一樣，不要輸入已經標准化的數據。順便也加載下edgeR這個R包
library(limma) library(edgeR) counts <- read.table(file = "conut_all.txt", sep = "\t", header = TRUE, row.names = 1, stringsAsFactors = FALSE) -
接着按照文檔的說明以及limma包的習慣，我們需要對count進行標准化以及轉化為log2的值，這里標准化的方法為TMM，使用edgeR里面的calcNormFactors函數即可
dge <- DGEList(counts = counts) dge <- calcNormFactors(dge) logCPM <- cpm(dge, log=TRUE, prior.count=3)這里prior.count值我粗略理解為是為了防止log2()遇到過於小的值
-
然后跟其他包一樣，設置分組信息
group_list <- factor(c(rep("control",2), rep("siSUZ12",2))) design <- model.matrix(~group_list) colnames(design) <- levels(group_list) rownames(design) <- colnames(counts) -
接下來就是常規的差異分析
fit <- lmFit(logCPM, design) fit <- eBayes(fit, trend=TRUE) output <- topTable(fit, coef=2,n=Inf) sum(output$adj.P.Val<0.05)
到這里為止，我們主要是用了limma包里對RNA-Seq差異分析的limma-trend方法，該方法主要適用於樣本間測序深度基本保持一致的情況，也就是說兩個樣本的文庫(reads數目)大小相差的不懸殊（說明文檔中是默認3倍以內？）
當文庫大小在樣本間變化幅度相當大的話，可以使用limma的voom方法，可按照下面的代碼流程：
-
count數據的輸入以及數據標准化還是跟之前的一樣
counts <- read.table(file = "conut_all.txt", sep = "\t", header = TRUE, row.names = 1, stringsAsFactors = FALSE) dge <- DGEList(counts = counts) dge <- calcNormFactors(dge) -
還是一樣需要分組信息
group_list <- factor(c(rep("control",2), rep("siSUZ12",2))) design <- model.matrix(~group_list) colnames(design) <- levels(group_list) rownames(design) <- colnames(counts) -
接下來進行voom轉化
v <- voom(dge, design, plot=TRUE)其實可以不進行TMM標准化，直接用count數據進行voom轉化，如：
v <- voom(counts, design, plot=TRUE) -
最后就是普通的差異分析過程
fit <- lmFit(v, design) fit <- eBayes(fit) output <- topTable(fit, coef=2,n=Inf) sum(output$adj.P.Val<0.05)
Summary
Limma長久以來就是一個非常流行的差異分析R包，其內容涉及的非常廣泛，用於RNA-Seq只是其內容的一小部分，並且使其處理RNA-Seq數據也使用芯片類似線性模型下，並且按照其說法，limma包比其他基於負二項式分布模型的差異分析R包更加的優秀。
其實差異分析不外乎數據的標准化以及統計模型分析差異兩個方面的作用，每個差異分析R包都有其自身的優點，個人理解，取舍在於自己的理解以及想法即可。
其實，自己對於limma包的理解還是比較粗淺的