Affymetrix芯片分析:获取差异表达基因系列二_Moderated t-test
library(limma)
library(tcltk)
library(affy)
filters <- matrix(c("CEL file", ".[Cc][Ee][Ll]", "All", ".*"), ncol = 2, byrow = T)
cellFiles <- tkSelectFiles(title = "Select CELs", multiple = TRUE, filters = filters, initialIndex = 1);
data.raw <- ReadAffy(filenames = cel.files)
design <- model.matrix(~-1+factor(c(1,1,2,2)))%构建实验设计矩阵(实验组2个和对照组2个)
colnames(design) <- c("control", "LPS") %实验组和对照组的列名称
eset.rma <- rma(data.raw)%调用RMA算法对数据进行预处理
eset<-exprs(data.raw)%构建样品的基因表达矩阵
fit <- lmFit(eset, design)%线性模型拟合
contrast.matrix <- makeContrasts(control-LPS, levels=design) %建立对比模型以分析两个实验条件下的基因表达差异
fit1<- contrasts.fit(fit, contrast.matrix)%根据对比模型进行差值计算
fit2 <- eBayes(fit1) %贝叶斯检验
results<-decideTests(fit2, method="global", adjust.method="BH", p.value=0.01, lfc=1.5)通过基于P-value的结果筛选获得差异表达基因
summary(results)
vennCounts(results)
vennDiagram(results)
heatDiagram(results,fit2$coef)
该代码段通过调用topTable函数对fit2进行分析,并设置相关参数以实现预期效果;具体设置包括:将系数设定为1;输出结果数量设为1万;采用Benjamini-Hochberg方法校正p值;按字母顺序对结果进行排序;并在必要时维持字母顺序
write.table(x, file="limma.xls", row.names=F, sep="\t")
y <- x[xP.Value < 0.01 & (xlogFC > 1.5 | xlogFC < -1.5& xAveExpr > 10),]
