GSVA分析

作者: 落寞的橙子 | 来源:发表于2020-04-03 03:03 被阅读0次

    GSVA其实就是pathway级别的差异分析,根据GSVA的原理其实就是计算每个sample的GSEA然后得出类似pathway enrich score,把这个可以当作芯片的表达数据一样,再用limma包分析差异基因。
    参考以下两个例子:
    使用GSVA方法计算某基因集在各个样本的表现
    充分理解GSVA和GSEA
    以及我一个例子:

    rm(list=ls())
    suppressMessages(library(GSVA))
    suppressMessages(library(GSVAdata))
    suppressMessages(library(GSEABase))
    suppressMessages(library(limma))
    
    #读入gmt文件,这个可以从MSigDB上下载,这边选的上gene symbol根据自己的data来选择
    gmt_file="~/c5.bp.v7.1.symbols.gmt"
    geneset <- getGmt(gmt_file)  
    
    #自行读入exp文件,我这边为行为gene symbol,列为样本(4个control,5个treatment样本)
    es <- gsva(as.matrix(exp), geneset,
                        min.sz=10, max.sz=500, verbose=TRUE)
    #得到gsva计算的数值后再用limma包做差异分析得到差异的pathway
    design <- model.matrix(~ factor(c(rep("cont",4),rep("treatment",5))))
    colnames(es) <- c("ALL", "contvstreatment")
    row.names(design)<-colnames(exp)
    fit <- lmFit(es, design)
    fit <- eBayes(fit)
    #这边是总的
    allGeneSets <- topTable(fit, coef="contvstreatment", number=Inf)
    #这边是差异的
    adjPvalueCutoff <- 0.001
    DEgeneSets <- topTable(fit, coef="MLLvsALL", number=Inf,
                           p.value=adjPvalueCutoff, adjust="BH")
    res <- decideTests(fit, p.value=adjPvalueCutoff)
    

    相关文章

      网友评论

        本文标题:GSVA分析

        本文链接:https://www.haomeiwen.com/subject/zjzqphtx.html