https://www.sohu.com/a/281843011_120055884
https://www.jianshu.com/p/cff8fe8b548a
之前介绍的一种是用R包做的图，利用的是log2foldchange 和ensemble ID
这次介绍的是GESA软件操作
原始分析数据
原始的FPKM值为
基因名（需要是ENTREZID格式）
构建其他物种
GSEA软件下载：https://www.gsea-msigdb.org/gsea/index.jsp 直接选择你操作系统对应的版本，下载安装即可。
准备上传文件
1.基因表达量表（.gct）
通常用.gct为后缀。文件第一行以“#1.2”开头；文件第二行的第一列为基因个数、第二列为样品个数；文件的第三行为表达谱的矩阵的title信息，第一列为基因symbol/探针号，第二列为基因／探针的描述信息，第三列以后为样品id。接下来的行对应每个基因／探针在每个样品中的表达信息。文件以tab作为分隔符。
image.png
第一列不能有重复
image.png
2. 样品表型分类文件(cls)——必需文件
样品表型分类文件需以.cls为后缀。文件第一行为三个数字，第一个是样品的总数，第二个是样品分为几类，第三个数字通常为1。第二行也通常三个字符串，第一个为#，第二个为分类1的名称，第三个位分类2的名称。第三行为每个样品的分类信息，0代表分类1，1则代表分类2。文件以空格或者tab分割。
GSEA 分析
导入数据
打开GSEA软件，点击左侧菜单栏的Load data选项，将数据按如下图所示方法导入
设置参数
按照下图标记的批注设置合适的参数。
Expression dataset（表达文件）： 选择上一步上传的表达gct文件
Gene sets database （功能基因集数据库）：GSEA包含了MSigDB数据库中的功能基因集，可以从中选择感兴趣的通路、癌症标记、转录因子数据库等。
Number of permutations（扰动／随机次数）：通常设置1000，此参数不可过小。
Phenotypes labels(样品表型分类文件)：选择上一步上传的表型cls文件
Collapse dataset to gene symbols：通常true
如果已经变为了gene symobols 则不需要进行转换
如果不是,则需要进行转换
Permutation type(扰动类型)： 通常选择phenotype，如果样品数目较少可以选择gene_set。
Chip platform(芯片类型)：如果表达gct文件的第一列为芯片探针id则此处需要选择对应的芯片平台，如果是基因symbol则无需选择。
image.png
image.png
image.png
gsea分析结束后
点击相应的success
跳到相应的分析结果
如果是其他物种,https://www.jianshu.com/p/cff8fe8b548a
可以构建心得基因集
<meta charset="utf-8">
基因集文件的准备相对比较复杂，我用到的是KEGG的通路基因集，也可以用GO的的基因集，甚至可以自己创造一个，只需要满足上图所示的文件格式。下面介绍一下如何从KEGG网站上下载某物种（小鼠）的基因集信息。
首先，在KEGG数据库中下载小鼠的KEGG通路数据库文件
大家可以看到，这个.keg文件和我们需要的.gmt文件还有不小的差异，那么，我们需要按如下方法，先在Linux下解析keg文件，然后用R制作.gmt文件：
#解析信息
perl -alne '{if(/^C/){/PATH:mmu(\d+)/;$kegg=$1}else{print "$kegg\t$F[1]" if /^D/ and $kegg;}}' mmu00001.keg > mmu00001_kegg2gene.txt
grep '^C' mmu00001.keg | grep "mmu" | sed 's/^C //g' | sed 's/\([0-9]*\)\s/\1\t/' > mmu00001_kegg2description.txt
#转换ID同时制作.gmt文件
#安装转换ID需要用到的各个R包
#安装data.table
install.packages("data.table")
#安装org.Xx.eg.db系列包中小鼠资源包
BiocManager::install("org.Mm.eg.db")
#安装clusterProfiler包，安装时需要“科学上网”
BiocManager::install("clusterProfiler")
#调用各个R包
library(data.table)
library(org.Mm.eg.db)
library(clusterProfiler)
#读入数据
path2gene_file<-"mmu00001_kegg2gene.txt"
tmp=read.table(path2gene_file,sep="\t",colClasses=c('character'))
dt <- tmp
dt$geneid <- rownames(dt)
# transform id
map_dt <- bitr(dt$V2, fromType = "ENTREZID",toType = c( "SYMBOL"),OrgDb = org.Mm.eg.db)
dt_merge <- merge(map_dt,dt, by.y = "V2", by.x = "ENTREZID")
# 可以把转换ID后的表输出
#write.table(dt_merge, file = "example42.txt", sep = "\t", col.names = T, row.names = T, quote = F)
# first column is kegg ID, second column is entrez ID
GeneID2kegg_list<- tapply(tmp[,1],as.factor(tmp[,2]),function(x) x)
kegg2symbol_list<- tapply(dt_merge$SYMBOL,as.factor(dt_merge$V1),function(x) x)
#输出数据，但描述列为NA
write.gmt <- function(geneSet=kegg2symbol_list,gmt_file='kegg2symbol.gmt'){
sink( gmt_file )
for (i in 1:length(geneSet)){
cat(names(geneSet)[i])
cat('\tNA\t')
cat(paste(geneSet[[i]],collapse = '\t'))
cat('\n')
}
sink()
}
write.gmt()
另外一种方法，构建斑马鱼物种的bp基因集
setwd("I:/ZHS_Seq_data/RUN")
library(msigdbr)
msigdbr_species() #可以看到，这个包涵盖了20个物种
zebrafish <- msigdbr(species = "Danio rerio")
zebrafish[1:5,1:5]
table(zebrafish$gs_cat) #查看目录，与MSigDB一样，包含9个数据集
#例如，本例中，我们要分析GO，因为mouse文件包含了所有的基因集，所以要查看GO在哪里，然后将需要的文件提出来。
table(zebrafish$gs_subcat)
zebrafish_GO_bp = msigdbr(species = "Danio rerio",
category = "C5", #GO在C5
subcategory = "GO:BP") %>%
dplyr::select(gs_name,gene_symbol)#这里可以选择gene symbol，也可以选择ID，根据自己数据需求来，主要为了方便
head(zebrafish_GO_bp)
zebrafish_GO_bp_Set = zebrafish_GO_bp %>% split(x = .$gene_symbol, f = .$gs_name)#后续gsva要求是list，所以将他转化为list
write.csv(zebrafish_GO_bp_Set,"zebrafish_GO_bp_Set.csv")
###names(zebrafish_GO_bp_Set)
##zebrafish_GO_bp_Set
write.gmt <- function(geneSet=zebrafish_GO_bp_Set,gmt_file='kegg2symbol.gmt'){
sink( gmt_file )
for (i in 1:length(geneSet)){
cat(names(geneSet)[i])
cat('\tNA\t')
cat(paste(geneSet[[i]],collapse = '\t'))
cat('\n')
}
sink()
}
write.gmt()
构好文件如下，第二列为NA
