ChIP-seq数据分析(一):从raw reads到peaks

ChIP-seq数据分析(一):从raw reads到peaks,第1张

最终的结果弊竖表明这6个数据的质量都很好,所以原文并没有进行额外的过滤 *** 作。需要注意的是,在该文献中ein2-5_air_MPGB_1和ein2-5_ethylene_MPGB_1两个阴性对照呈现出reads重复租携大率高的特征

-M 1 --best --strata 的作用是将多比对中得分最高的一条比对记录保留,原文在这里采用的是保留唯一比对的reads记录,有些区别。

bowtie比对完之后会在log文件中报告比对率,比如

还可以用samtools flagstat简单统计比对率

又统计了一下fq里面的reads数量

这一步没有什么实质性的作用,但能加深对reads比对率、唯一比对、多比对、reads数、比对记录数等的理解。

"Duplicate reads were removed using SAMtools."相关的NC文章是这样做的。

并用samtools index为4个待分析的bam建立索引。

文章只强调了--nomodel, -p 0.01这两个参数。

看参数的值就能知道bw存储的是什么信息:横坐标是在基因组上的一对起始位置,窗口大小是50bp,纵坐标是将深度标准化之后得到的RPKM。

除了bamCoverage,bamCompare也能将bam->bw,并且同时考虑处理和对照,以消除噪声。原文是这样说的: To show ChIP binding signal surrounding TSSs or in gene bodies, read coverage was first calculated using the bamCompare tool (RPKM, Log2(ChIP/IgG) in deepTools.

在igv下面看看bed, bam, bw, bgd是什么样子

可以看到,bam, bw, bdg整体趋势一致,并且跟对照相比,可以显示出peaks。

针对的是summits.bed文件,看peaks(一个峰,精确到了单碱基位置)属于哪一个基因?或是离哪一个基因最近?属于哪一个区域?

这个summits的注释跟变异检测注释原理差不多,只隐陵是分区不一样。

前面的帖子 一文阐述链特异性测序——stranded? reverse-stranded? un-stranded? 分享了链特异性测序的相关内容,今天继续分享一个内容: 如何判断测序数据是否是链特异性?

查找原文描述,一般在 method 中,由此判断瞎者态是否是链特异性测序。

比对完成后,将 bam 文件在 igv 上进行可视化,以单端测序为例,如果你发现测序的read方向总是与所在基因的方向一致或者总是相反,那么就是链特异性测序;反之,则不是。

利用 RseQC 进行判断,其python小脚本 infer_experiment.py 可以帮助我们判断链特异性测序。

安装完成后发现,其实并没有安装软件,而是多了一些脚本,其中有一个就叫 infer_experiment.py

部分参数:

(1) -i 比对生成的bam文件(可以不用排序)

(2) -r gtf转bed12文件产生的bed文件。 技能——gtf转为bed12

(3) -s 从所有的reads中抽取多少进行统计(默认200k)

(4) -q unique map的mapq阈值

你会发现有接近80%的reads都是比对到哪条链,基因就在哪条链,占绝对优势。 这个时候我们有理由认为这个单端测序的数据就是链特异性测序的数据。

详细来看: ++ 中第一个 + ,表示这个read比对上正链,第二个 + ,表示这个read所在的基因也是位于正链的。

这个双端测序的数据显然是磨源非链特异性的,详细来看:

1++,1--,2+-,2-+ 表示read1比对上的链和基因所在的链一样,read2比对上的链和基因所在的链是相反的;

1+-,1-+,2++,2-- 表示read1比对上的链和基因所在的链是相反的,read2比对上的链和基因所在的链一样。嫌隐

问题就在于这两种情况的占比是几乎一样的,说明read的链与基因所在的链完全是随机的,显然就不是链特异性测序了。

完结,撒花~

linux学的半半拉拉,觉得需要用实战来激励一下,之后再回过头还要再学linux的哦。所以先来ChIP-seq吧,这个挺急需的,去年的数据估计也需要用这一套流程来跑。感谢Jimmy大神无私地录教程。 https://www.bilibili.com/video/BV16s411T7Fh?spm_id_from=333.999.0.0

raw data转换,quality control。

call peak-通常是用MACS2来找到peaks。有的可以IGV可视化。

包括结合峰相关基因的注释、GO分析KEGG通路富集分析和motif分析,主要是 motif和peaks的注释。

包括差异结合峰相关基因的注释、GO分析、KEGG通路富集分析和motif分析。即多次重复的相关性,或者处理前后peaks差异的对比。

还能对目的蛋白的调控机制做进一步的深吵搭入分升局拿析。

从知乎大神-生信宝典这里学习到 https://zhuanlan.zhihu.com/p/436823899

这些之前都知道啦。bam和sam文件可以帮助我们探索reads在参考基因组中的比对情况。

只包含区域和区域的覆盖度信息,文件更小,用于可视化更方便,可以导入基因组浏览器(Genome Browser)中进行可视化,以查看reads在参考基因组各个区域的覆盖度并检测测序深度。这几个文件在[ChIP-seq数据分析]、Call Peak阶段会生成,可以利用[IGV]等工具进行可视化,方便查看组蛋白修饰的程度。

wig文件包括染色体长度,步长多少,span多少

bw文件,bigWig是 wig文件的二进制压缩格式,可通过wigToBigWig工具转换。bw文件可以直接导入IGV。

bedGraph: bed与wig的结合,更省空间和灵活,展示信息与wig类似。

分析过程中的bed文件一般代表区域信息,如表示Peak位置的bed文件,表示基因注释的bed12文腊旅件。

下文可以详读。

Pedro Rosmaninho et al, Zeb1 potentiates genome‐wide gene transcription with Lef1 to promote glioblastoma cell invasion, The EMBO Journal (2018). DOI: 10.15252/embj.201797115


欢迎分享,转载请注明来源:内存溢出

原文地址: http://outofmemory.cn/tougao/8227372.html

(0)
打赏 微信扫一扫 微信扫一扫 支付宝扫一扫 支付宝扫一扫
上一篇 2023-04-14
下一篇 2023-04-14

发表评论

登录后才能评论

评论列表(0条)

保存