首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >引用2000多次的ATAC经典文献也在用的peak calling软件-Genrich

引用2000多次的ATAC经典文献也在用的peak calling软件-Genrich

作者头像
生信修炼手册
发布2020-05-07 16:14:48
9800
发布2020-05-07 16:14:48
举报
文章被收录于专栏:生信修炼手册生信修炼手册

在之前的文章中,我们解读过一篇引用达2000多次的ATAC经典文献

引用2115次的ATAC经典论文解读

这篇文章中,使用了Genrich这个软件来进行peak calling。该软件适用于chip_seq, DNase_seq, ATAC_seq等多种文库的peak calling,源代码保存在github上,链接如下

https://github.com/jsh58/Genrich

该软件用法简单,只需要输入排序之后的bam文件即可,其peak calling的基本算法如下所示

  1. 对原始bam文件进行过滤,该软件只针对双端测序的bam文件,默认情况下只保留双端同时比对到基因组上的序列。和macs2不同,该软件支持multimapping reads,对于比对到基因组多个位置的reads, 每个位置会计算一个权重。同时还可以过滤PCR重复序列,去除指定染色体上的序列
  2. 如果提供了control样本,首先利用control样本构建一个背景
  3. 对于treat样本,利用log-normal分布来计算每个位点的p值
  4. 多重假设检验的校正,将p值转换为q值
  5. 计算统计学显著,即p值或者q值小于0.05的区域,对应的曲线下的面积,即AUC, 设定AUC的阈值,如果一个区域的AUC大于阈值,则定义为peak

该软件有两种运行模式,第一种适用于chip_seq, 第二种适用于ATAC。对于ATAC的reads,需要对reads比对位置进行shift, 所以peak calling的模式也会有所不同,图示如下

默认情况下为第一种模式,通过-j参数可以调整为ATAC模式。默认模式下,将fragment两个末端对应的区域直接作为peak, 以chip_seq为例,抗体抓下来的reads就是蛋白结合区域的reads。而ATAC文库中Tn5切割的序列是在结合区域的两侧,所以在ATAC模式中,以fragment的中心点为基准,来判断真实的peak区域。

对于ATAC数据的peak calling而言,该软件的用法如下

Genrich \
-t SRR5427886.bam \
-o SRR5427886.narrowPeak  \
-f SRR5427886.log \
-r \
-x \
-e chrM \
-E wgEncodeDukeMapabilityRegionsExcludable.bed.gz \
-q 0.05 \
-a 20.0

-t参数指定输入的bam文件,-o参数指定输出的peak文件,-f参数指定log文件的名称。-r参数表示去除PCR重复,-x参数表示保留只有一端比对上基因组的reads,-e参数剔除指定染色体上的序列,-E参数剔除bed文件中指定的基因组区域。-q参数指定qvalue的阈值,-a参数指定AUC面积的阈值,

对于ATAC的peak calling而言,除了MACS2, 该软件也是一个不错的选择。

本文参与 腾讯云自媒体分享计划,分享自微信公众号。
原始发表:2020-01-11,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 生信修炼手册 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体分享计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档