生物信息学入门 GEO芯片数据差异表达分析时需要log2处理的原因

首先借用一张图,通常使用limma处理时,需要经过log2后的矩阵作为表达矩阵输入。根据log2FC的定义,这个数字表示变化倍数经过log2后的一个值,比如log2FC=1,则变化为2倍;log2FC=2,则变化为4倍。这是常用的一种表述方法。

在使用limma函数计算时,如果输入的矩阵没有经过log2处理,则会把FC当成log2FC输入,这或许是因为limma默认输入的是log2后的表达式。这里有必要提到log的一个运算,即,

可见对于已经log2后的数据,计算log2FC = log2(A/B)只需要直接使用log2A-log2B。所以如果给出的是一个未经log2的数值,函数也会直接相减以得到log2FC,这就导致计算出来的差异表达高达几百甚至上千。

并且,通过RMA法进行预处理时,已经经过了log2。

但是在GEO中,一些Series Matrix File(s)仍是没有log2或者标准化,关于判断方法见下帖:

GEO芯片数据差异表达分析时是否需要log2以及标准化的问题

https://blog.csdn.net/tuanzide5233/article/details/88542558

使用GEO芯片数据通过limma包进行差异表达的教程

https://blog.csdn.net/tuanzide5233/article/details/83541443

 

 

差异表达分析是一种常用的方法,用于比较不同条件下基因表达水平的变化。GEOquery和limma是R语言中广泛使用的两个包,可用于处理分析芯片数据差异表达。下面是使用GEOquery和limma进行差异表达分析的步骤: 1. 下载和导入芯片数据 使用GEOquery包中的getGEO函数下载芯片数据并导入到R中。例如,如果您要下载GSE12345数据集,可以使用以下代码: ``` library(GEOquery) gset <- getGEO("GSE12345") ``` 2. 数据质量控制 在进行差异表达分析之前,需要数据进行质量控制。使用GEOquery包中的summary函数和plotPCA函数可以对芯片数据进行基本的质量控制。例如,可以使用以下代码绘制PCA图: ``` library(limma) library(GEOquery) gset <- getGEO("GSE12345") edata <- exprs(gset[[1]]) edata <- t(edata) edata <- na.omit(edata) fit <- prcomp(edata) plotPCA(fit) ``` 3. 数据处理芯片数据进行归一化和标准化处理,以消除不同芯片之间的差异,并确保数据符合正态分布。常用的预处理方法包括RMA、GCRMA、MAS5等。使用limma包中的normalizeBetweenArrays函数可以对芯片数据进行预处理。例如,可以使用以下代码对芯片数据进行RMA预处理: ``` library(limma) library(GEOquery) gset <- getGEO("GSE12345") edata <- exprs(gset[[1]]) edata <- t(edata) edata <- na.omit(edata) edata <- backgroundCorrect.RMA(edata) edata <- normalize.quantiles.RMA(edata) edata <- log2(edata) ``` 4. 差异表达分析 使用limma包中的lmFit函数和eBayes函数可以进行差异表达分析。lmFit函数用于拟合线性模型,eBayes函数用于对差异表达结果进行统计显著性检验。例如,可以使用以下代码进行差异表达分析: ``` library(limma) library(GEOquery) gset <- getGEO("GSE12345") edata <- exprs(gset[[1]]) edata <- t(edata) edata <- na.omit(edata) edata <- backgroundCorrect.RMA(edata) edata <- normalize.quantiles.RMA(edata) edata <- log2(edata) factors <- c(0,0,0,1,1,1) design <- model.matrix(~factors) fit <- lmFit(edata, design) fit <- eBayes(fit) results <- topTable(fit, adjust="BH", sort.by="P", n=1000) ``` 上述代码中,factors表示芯片数据中不同样本的分组信息,design表示设计矩阵,fit表示拟合的线性模型,results表示差异表达结果。 5. 结果分析 根据差异表达结果,可以进行进一步的功能分析、通路分析等。常用的工具包括ClusterProfiler、GOstats、KEGGprofile等。例如,可以使用以下代码进行GO分析: ``` library(ClusterProfiler) gset <- getGEO("GSE12345") edata <- exprs(gset[[1]]) edata <- t(edata) edata <- na.omit(edata) edata <- backgroundCorrect.RMA(edata) edata <- normalize.quantiles.RMA(edata) edata <- log2(edata) factors <- c(0,0,0,1,1,1) design <- model.matrix(~factors) fit <- lmFit(edata, design) fit <- eBayes(fit) results <- topTable(fit, adjust="BH", sort.by="P", n=1000) genes <- rownames(results) geneList <- names(genes) geneList <- names(genes)[abs(results$logFC) > 1 & results$adj.P.Val < 0.05] ego <- enrichGO(geneList, OrgDb="org.Hs.eg.db", ont="BP") barplot(ego) ``` 上述代码中,geneList表示差异表达基因列表,ego表示进行GO分析所得到的富集结果。
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值