使用ChIPSeeker进行ChIP-seq, ATAC-seq,cut&tag等富集峰的基因组注释

二代测序产生的数据类型

常规的下一代高通量测序(next generation sequencing, NGS)实验通常产生大量短片段(reads),通常我们需要将这些reads比对到参考基因组/转录组上,即将它们置于生物学上有意义的基因背景下,才能获得有意义的结果。一般我们认为会产生两种类型的数据(当然两者并无严格意义上的区分):

1.表达类

一般为固定区域,关注于定量比较。例如转录组测序结果中,数据库中已有mRNA基因的表达,lncRNA基因的表达等,这类结果一般以矩阵形式存储,第一列是名字,其余列是表达值。CircRNA表达,miRNA表达等归于此类。

2.基因组区域(富集峰)

一般区域不固定,关注于定性。例如ChIP-seq、ATAC-seq、Cut&tag等比对后获得的富集峰。一般以bed格式存储,第一列是染色体,第二列是富集峰起始坐标,第三列是富集峰终止坐标(图1)。eccDNA,m6A,MeDIP等归于此类。

图1. 表达类(区域固定) vs 富集峰类(区域不固定)

什么是富集峰注释?

基于抗体富集的原理,众多reads片段比对到基因组上某区段,会形成一个类似山峰的富集区。由于我们是在基因组背景下进行生物医学研究的,因此需要将基因组区域(富集峰,peak)与基因联系起来,即确定峰落在哪个基因上,落在该基因的哪种基因组特征上,距离TSS的位置是多少bp等,然后才能进行后续的功能研究。这个过程叫做富集峰注释(peakannotation,图2)。如果仅关注某几个区域,就不需要用软件注释,建议直接用IGV或者UCSC Genome Browser查看。

富集峰注释的难点

虽然CHIP-Seq已经有近20年的历史,然而由于基因组上基因的结构非常复杂,不同注释软件在具体细节处理上往往不同,从而导致同一批数据,用不同的软件进行注释,获得的结果略有不同(大同小异)。

1.注释到转录本还是注释到基因?

由于一个基因可能包含多个转录本,因此,我们在注释的时候,到底是注释到基因水平还是注释到转录本水平?

2.基因位置重叠怎么处理?

由于同一个位置可能存在多个基因,如果两个基因的坐标有重叠,我们到底是注释到A基因还是注释到B基因?

3.最邻近怎么判断?

如果富集峰的中点(或者顶点)正好落在两个基因的中间,那么这个富集峰是注释到A基因还是注释到B基因,如何定义最邻近?

4.基因组特征如何分类?

不同文章、软件对基因组特征分类不同,例如有的分为promoter、intron、exon,5’UTR和3’UTR;有的分为:upstream、promoter、intron、exon、downstream等

5.如何定义?

不同文章、软件定义promoter区也不同,有的定义TSS上游3K到下游3K都是promoter,有的定义TSS上游200bp到下游800bp是启动子。

6.不同基因组特征的注释优先级

当一个很长的富集峰横跨同一个基因的intron、exon、3’UTR时候,这个富集峰该分到什么特征中呢?

7.注释数据库版本问题

同样是human hg38,如果注释库版本不同,那么注释结果也会有差异。原因是:虽然基因组序列不变,然而注释库却更新频繁,有基因会更新坐标,有基因会添加新的转录本,有基因会从非编码基因变成编码基因等。

图2. 富集峰注释及难点

凡此种种,给我们的注释工作带来了巨大困难。然而,作为用户(调包侠),我们基本不用深究注释背后的细节。我们需要做的就是:找一个引用比较多的注释工具,默认参数进行注释即可。

常见富集峰注释软件

表1. 常见富集峰注释工具

为什么要用新版注释?

由于注释数据库频繁更新,如果你使用的注释还是N年前的,那么reviewer在公共数据库(例如UCSC、Ensembl、NCBI)上使用网站默认版本查询时,就有可能查不到你的基因,或者你N年前的数据,与新的数据联合分析时,由于使用的注释数据库不同,取交集时,会漏掉一些基因。因此,我们强烈建议所有的测序数据,包括RNA-seq、ChIP-seq、m6A-seq等都使用同一套注释库进行注释分析,并在结果中明确说明所使用的注释库版本。这对于在不同公司,不同时间做的测序结果来说,是非常重要的。

由于上述所列在线工具都是N年前的,所以我们使用ChIPSeeker R包搭建了一个简易的在线peak注释工具,可以对人、大鼠、小鼠的ChIP-seq,ATAC-seq,cut&tag等富集峰进行一键注释。

1,打开绘图页面

首先,使用浏览器(推荐chrome或者edge)打开ChIP-Seq富集峰注释页面。左侧为常见作图导航,中间为数据输入框和可选参数,右侧为描述和结果示例。也可以在搜索框中搜索peak,找到注释页面。

http://www.bioinformatics.com.cn/basic_chipseq_atacseq_peak_annotation_by_chipseeker_t017

图3. 注释页面

2,示例数据

点击右侧“示例数据”链接下载excel格式的示例数据。

示例数据包括4列,分别为:chr,start,end,-LOG10(pvalue)。

注意:为了遵循各大数据库的使用,这里染色体必需使用chr+数字,即:chr1-22、chrX、chrY、chrM等。

图4. 输入数据示例

3,粘贴示例数据

直接拷贝示例数据中的ABCD四列数据,然后粘贴到输入框。注意必需带每一列的说明行(header),此行将用于最终的excel表头。

注意:不是拷贝excel文件,是拷贝excel文件里边的数据。另外粘贴到输入框后,格式乱了没关系,只要在excel中是整齐的就行。同时数据矩阵中不能有空的单元格,中文字符等。

图5. 必需输入

4,修改参数,并提交

我们设置了promoter区的范围选项,及注释所用的物种及注释版本选项(例如human所使用的是hg38基因组,注释版本为Ensembl v108),当前仅支持human,mouse和rat。后续将支持更多物种。

图6. 可选参数

5,提交出图

粘贴好输入数据,调整好参数(重点是物种及注释版本)后,点击提交按钮,约30秒钟后(取决于数据多少),会在页面右侧出现peak分类饼图及excel格式的数据下载链接,请下载后解压查看。

图7. 结果页面


图8. 注释结果

结果说明

以peak_class为界,结果包括两部分:左侧为输入的内容,其中start添加了1 bp(因为bed格式是0-based,这里变成了1-based),并添加了peak长度信息(end-start+1);右侧为注释信息,包括:peak分类,基因位置,基因/转录本注释等信息。并基于peak_class的数据绘制了peak分布饼图。

注意:

1,由于peak注释与注释库及优先级关系密切,因此最终放在paper里边的图,以IGV可视化结果为准。

2,输入peak默认不考虑链,如需更精细地注释,请参考ChIPSeeker R包。

3,默认一个peak仅注释到一个转录本

没有预览就是没有出图/结果,这时请参考示例数据,检查输入数据的格式。或者使用我们提供的小工具pyinstaller打包python脚本为exe可执行文件实例:错误排查小脚本检查输入。

参考文献:

1. Yu G, Wang LG, He QY. ChIPseeker: anR/Bioconductor package for ChIP peak annotation, comparison and visualization.Bioinformatics. 2015 Jul 15;31(14):2382-3. doi: 10.1093/bioinformatics/btv145.Epub 2015 Mar 11. PMID: 25765347.

2. Huang W, Loganantharaj R, Schroeder B, FargoD, Li L. PAVIS: a tool for Peak Annotation and Visualization. Bioinformatics.2013 Dec 1;29(23):3097-9. doi: 10.1093/bioinformatics/btt520. Epub 2013 Sep 4.PMID: 24008416; PMCID: PMC3834791.

3. Heinz S, Benner C, Spann N, Bertolino E, LinYC, Laslo P, Cheng JX, Murre C, Singh H, Glass CK. Simple combinations oflineage-determining transcription factors prime cis-regulatory elementsrequired for macrophage and B cell identities. Mol Cell. 2010 May28;38(4):576-89. doi: 10.1016/j.molcel.2010.05.004. PMID: 20513432; PMCID:PMC2898526.

4. McLean CY, Bristor D, Hiller M, Clarke SL,Schaar BT, Lowe CB, Wenger AM, Bejerano G. GREAT improves functionalinterpretation of cis-regulatory regions. Nat Biotechnol. 2010May;28(5):495-501. doi: 10.1038/nbt.1630. Epub 2010 May 2. PMID: 20436461;PMCID: PMC4840234.

5. Zhu LJ, Gazin C, Lawson ND, Pagès H, Lin SM,Lapointe DS, Green MR. ChIPpeakAnno: a Bioconductor package to annotateChIP-seq and ChIP-chip data. BMC Bioinformatics. 2010 May 11;11:237. doi:10.1186/1471-2105-11-237. PMID: 20459804; PMCID: PMC3098059.

微生信助力高分文章,用户73000+,引用990+

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 194,088评论 5 459
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 81,715评论 2 371
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 141,361评论 0 319
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 52,099评论 1 263
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 60,987评论 4 355
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 46,063评论 1 272
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 36,486评论 3 381
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 35,175评论 0 253
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 39,440评论 1 290
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 34,518评论 2 309
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 36,305评论 1 326
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 32,190评论 3 312
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 37,550评论 3 298
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 28,880评论 0 17
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,152评论 1 250
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 41,451评论 2 341
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 40,637评论 2 335

推荐阅读更多精彩内容