0055-【生物数据库】-如何进行RNA差异基因KEGG注释分析-kobas在线分析

1. 有参物种使用gene ID的方法

1. 差异基因文件准备

只需要用到两列

• ENTREZ_GENE_ID
• logFC

geneNames ENTREZ_GENE_ID normalAve tumorAve logFC pValue qValue CCL23 6368 95.0 5. -4.0 2.07E-31 5.99E-29 COLEC10 10584 1459. 83. -4. 2.11E-31 6.00E-29 FAM189B 10712 383. 1289. 1. 2.17E-31 6.08E-29 CDC45 8318 12. 258. 4. 3.59E-31 9.94E-29 RCAN1 1827 11046.97758 2309. -2. 3.90E-31 1.07E-28 N4BP2L1 90634 2644.65753 734.73331 -1. 4.57E-31 1.23E-28 FCN3 8547 6777. 389. -4. 5.41E-31 1.44E-28 UHRF1 29128 15. 327. 4. 5.73E-31 1.50E-28 HMMR 3161 25. 407. 4.00 8.18E-31 2.12E-28 NEK2 4751 18. 390. 4. 9.48E-31 2.43E-28

选择基因的ID作为输入文件

6368 10584 10712 8318 1827 90634 8547 29128 3161 4751

2. 登陆kobas数据库

网站：http://kobas.cbi.pku.edu.cn/

进入 Gene-list-Enrichment
http://kobas.cbi.pku.edu.cn/anno_iden.php

输入数据类型：

• Fasta Protein Sequence ——蛋白序列
• Fasta Nucleotide Sequence——核酸序列
• Tabular BLAST Output——blast输出的表格
• Entrez Gene ID——基因ID
• UniProtKB AC
• Refseq Protein ID
• Ensembl Gene ID

3. 选择

1. 输入类型选择：Gene ID

2. 物种选择：Homo sapiens (human)

3. 粘贴Gene ID列表

4. 数据库 Clear All取消Pathway、Disease、GO全部选项，只选择KEGG Patway

点击RUN

4. 在线分析完成，输出结果

5. 输出文件说明

统计学检验方法：超几何检验、FIsher精确检验
FDR校正方法：Benjamini and Hochberg，需要补充此方法

Statistical test method: hypergeometric test / Fisher's exact test FDR correction method: Benjamini and Hochberg

输出表格：

• Term KEGG的注释类
• Database 数据库类型
• ID Term对应的ID
• input number 富集到这个Term的输入基因个数
• Background number 数据库中富集到这个通路的总有基因数量
• P-value P值
• Corrected P-Value 校正后P值
• Input 输入的Gene ID，如果多个，以|号分开
• Hyperlink 网页链接

图片会将对应的Gene name标志为红色

6. 软件安装准备

由于bioconductor外网链接慢，使用conda的方法安装，同时安装依赖的包

conda install bioconductor-clusterprofiler 

7. 画图

# 初始化环境 rm(list=ls()) # 安装软件 #source("https://bioconductor.org/biocLite.R") #biocLite() #biocLite("clusterProfiler") #biocLite("pathview") # 设置通路 setwd("/home/toucan/Project/001.kegg_map") # 加载库 library("clusterProfiler") # 读入文件，不检测name rt=read.table("input.txt",sep="\t",header=T,check.names=F) rt # 构建gene id为行名称的，logFC  geneFC=rt$logFC geneFC gene <- rt$ENTREZ_GENE_ID gene names(geneFC)=gene geneFC #kegg # 保存输出文件 # 设定物种，qvalue小于0.05才输出，readable是否输出转换为gene name kk <- enrichKEGG(gene = gene, organism = "human", pvalueCutoff = 0.05,qvalueCutoff = 0.05) class(kk) kk as.data.frame(kk) write.table(as.data.frame(kk),file="KEGG.xls",sep="\t",quote=F,row.names = F) # 生成barplot pdf(file="KEGG.barplot.pdf") barplot(kk, drop = TRUE, showCategory = 12) pdf(file="KEGG.cnetplot.pdf") # 生成网络图，需要通路描述列、输入gene ID列组成 #cnetplot(kk,categorySize = "geneNum", foldChange = geneFC) library("pathview") keggxls=read.table("KEGG.xls",sep="\t",header=T) # 联网，将map图片下载 for(i in keggxls$ID){ pv.out <- pathview(gene.data = geneFC, pathway.id = i, species = "hsa", out.suffix = "pathview")}

7. 输出结果

输出富集的表格：

ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count hsa04110 Cell cycle 19/199 124/7431 5.54E-10 1.37E-07 1.31E-07 8318/7272/890/1870/701/4085/4998/4171/4175/898/23594/1031/4172/4616/8317/4176/4174/9134/993 19 hsa03030 DNA replication 11/199 36/7431 1.29E-09 1.60E-07 1.53E-07 2237/4171/4175/10535/5984/4172/5558/5424/23649/4176/4174 11 hsa03440 Homologous recombination 8/199 41/7431 1.03E-05 0.000 0.000 /8438/5888/7517/5424/641/7516/25788 8 hsa05222 Small cell lung cancer 11/199 93/7431 3.44E-05 0.00 0.00 1870/898/3910/4616/1282/3655/1284/9134/5743/3915/1163 11 

每个通路生成三个文件

• hsa03030.pathview.png
• hsa03030.png
• hsa03030.xml

输出富集的图片

输出伏击通路下载的map：
有差异基因显示，红色为正相关，绿色为负

同时，输出网站原始下载的，无颜色标注

非模式生物

以序列作为输入文件

>seq1 CTAATTTTGATGTAACAATAAGCAAATCCATCTCATTGACATGTCAACTTACCTTAATCTTTAATAAGTG ATAAAGTCATATGTATGCCAAAAATTGCCTTAGCATTGCGTTATGACCTACCGTTAGTAGATGTCTGATT >seq2 AGTCTCGAATACAACTTGTTGCTGCGCGGACGCGAATCGCTCAGTACGGACGTCTTGAGCTCGAATCCTC GGCCATATCTGTGCTCTCGATCGCAGCGTTTGCTAATTCGAAGATCGTGCTAATCGAAGTACCGAGAAAT 

注意，物种应选择KO,但会笔记慢

If you choose KO, Please input no more than 200 lines at one time.

运行中：

http://kobas.cbi.pku.edu.cn/wait_kobas.php?taskid=9220 Your task is still running, your task id is 9220, you can get the results automatically when the task is finished. Also you can use the task id to fetch results at the result retrive page in the future.

等待输出