干货 | GSVA（Gene Set Variation Analysis）基因集变异分析

2022-08-20 06:47:20, 欧易生物 上海欧易生物医学科技有限公司

GSVA基因集变异分析，是在无参数且无监督的条件下，将Microarray和RNA-seq数据基因集进行富集的分析方法。GSVA能将一个基因-样本的数据矩阵（microarray data, FPKM, RPKM等）转换成基因集-样本矩阵。基于该矩阵可以进一步分析各个样本中基因集（如KEGG通路）的富集情况。由于GSVA的结果为一个基因集-样本的富集矩阵，相比其它基因集富集方法如GSEA （Gene Set Enrichment Analysis），下游分析会有更大的自由度。

在下游分析中，常见的分析是基于R包limma的差异分析。通过差异分析，找出实验组与对照组间存在富集值差异的通路。

准备：

R（推荐在RStudio上操作）、R包GSVA、GSEABase和limma


输入文件：

表达矩阵文件（FPKM、counts、RPKM等）和gmt矩阵文件（GSEA分析输入的通路矩阵文件）


表达矩阵文件示例：


gmt矩阵文件示例：




代码示例

# 依赖包安装

install.packages("BiocManager")BiocManager::install("GSVA")library(''GSVA'')BiocManager::install("GSEABase")library(''GSEABase''


#输入fpkm.xls 表达矩阵文件

expset  <- ''fpkm.xls''expset <- read.table(expset, sep = ''\\t'', header = TRUE, row.names = 1)


#Gaussian方法需要对表达矩阵log化

expset <- log2(expset + 1)matrix_expset <- data.matrix(expset)


# keg.backgroud.gmt gmt矩阵文件输入

gmt <- ''kegg.backgroud.gmt''geneSet <- getGmt(gmt)


# method可调整三种算法，针对不同需求，默认为gsva

kcdf参数包含不同统计方法，针对不同类型的输入矩阵应采用不同方法。与之同时纳入分析最小基因集，最大基因集也是可调的。输出结果将基因-样本的表达矩阵转换为通路-样本的富集分数矩阵。

gsva_scores <- gsva(matrix_expset, geneSet, method = ''gsva'', min.sz = 15 , max.sz = 5000, kcdf = ''Gaussian'', parallel.sz = 10)


# 接着利用limma进行差异分析

BiocManager::install("limma")library(''limma'')groups <- factor(c(rep(''A'', 3,), rep(''B'',3)))design <- model.matrix(~ groups + 0)rownames(design) <- colnames(gsva_scores)compareE <- makeContrasts(groupsA-groupsB, levels = design)fit <- lmFit(gsva_scores, design)fit <- contrasts.fit(fit, compareE)fit <- eBayes(fit)


#得到未筛选的limma两组间比较结果

fit_all <- topTable(fit, coef = 1, adjust = ''BH'', n = Inf)fit_all <- cbind(rownames(fit_all),fit_all)colnames(fit_all) = c("geneset", "logFC", "avgExp", "t", "pval", "FDR", "B")


#于是便得到了一个gsva富集矩阵及组间的通路差异结果。


最后，可以对差异结果做一定的可视化，如热图或者柱状图。


#筛选部分结果

diff_result = fit_all[fit_all$pval < 0.05,]diff_result$up_down = ifelse(diff_result$t<0,''Down'',''Up'')


#做柱状图

library(ggplot2)pp = ggplot(diff_result, aes(reorder(geneset, t), t))  + geom_col(aes(fill=up_down)) + scale_fill_manual(values=c("#6CC570","#2A5078")) +  coord_flip() +    labs(x="Gene set", y="t value of GSVA Score" ) +    theme_minimal() +    geom_text(data = diff_result[diff_result$up_down == ''Down'',], aes(x= geneset, y=0, label = geneset), hjust = 0, size = 2.5)+    geom_text(data = diff_result[diff_result$up_down == ''Up'',], aes(x= geneset, y=0, label = geneset), hjust = 1, size = 2.5)+ theme_bw() +    theme(axis.line.x = element_line(colour = "black"))+    theme(axis.line = element_blank()) +    theme(axis.text.y=element_blank()) +theme(panel.grid =element_blank(),panel.border = element_blank(),axis.text.x=element_text(size=6),legend.text=element_text(size=8),axis.ticks = element_blank())+ labs( fill = "Up_Down")


#存图

ggsave("barplot.pdf",         height = 10,width = 8,plot = pp, limitsize = FALSE)


- OE BIOTECH

OE BIOTECH -



END

排版人：小久


原创声明：本文由欧易生物(OEBIOTECH)学术团队报道，本文著作权归文章作者所有。欢迎个人转发及分享，未经作者的允许禁止转载。

点击“阅读全文” 收获更多精彩

  • 客服电话： 400-6699-117 转 1000
  • 京ICP备07018254号
  • 电信与信息服务业务经营许可证：京ICP证110310号
  • 京公网安备1101085018
  • 客服电话： 400-6699-117 转 1000
  • 京ICP备07018254号
  • 电信与信息服务业务经营许可证：京ICP证110310号
  • 京公网安备1101085018

Copyright ©2007-2024 ANTPEDIA, All Rights Reserved