ChIP-seq | ATAC-seq | Bulk | 新手指南 | 重复的处理

 

2023年05月17日

cut&run call出了太多的peak，大概6w个，但TF的丰度不如Histone mark，所以这些peak的质量是会被人质疑的。虽然我是把replicate合在一起call的，理论上的peak是更靠谱的，但到IGV里可视化时效果就没那么好了，除非你把bigwig文件合并一下，但审稿人还是想看你的replicate的一致性情况，作为一个重要的QC指标。

如此看来ChIP-seq还是有其不可取代的地位的，因为TF的丰度问题，Histone则肯定用cut&run更好，噪音更少，还是个信噪比问题。

关于replicate重复的处理

如果只有两个重复样本，那就用bedtools的intersect功能，最好输出一个peak的完整peak region，不要截断；

如果有三个及以上的重复样本，那就用encode的IDR工具，即一个peak有超过两个样本覆盖就算高质量的peak；

# intersect
# sort -k8,8nr # sort by score
bedtools intersect -wa -a HDAC1-1.sorted.mapped.bam_peaks.narrowPeak -b HDAC1-2.sorted.mapped.bam_peaks.narrowPeak | grep chr | sort -k1,1 -k2,2n > HDAC1.narrowPeak
bedtools intersect -wa -a HDAC2-1.sorted.mapped.bam_peaks.narrowPeak -b HDAC2-2.sorted.mapped.bam_peaks.narrowPeak | grep chr | sort -k1,1 -k2,2n > HDAC2.narrowPeak

# merge
cat HDAC1.narrowPeak HDAC2.narrowPeak | cut -f1-3 | sort -k1,1 -k2,2n | bedtools merge -i - > HDAC1_2.merge.narrowPeak

 

参考：

• Handling replicates with IDR - HBC
• Combine ChIP-seq peaks from multiple replicates via consensus voting
• bedtools intersect

 

---

 

没想到我是先玩了Cut&Run和单细胞ATAC-seq，搞通了再来分析ChIP-seq | ATAC-seq，因为之前接触到的数据太烂了，所以什么有意义的东西都没搞出来，又没有重复，导致分析和评估无从入手。

现在这个lab还是比较正规的，数据有小瑕疵，但整体而言质量是非常靠谱的，各种数据之间也能完美的自我验证，是非常适合初学者来入门的。

 

Cut&Run的分析我基本吃透了，因为没什么背景噪音，IgG用spike-in校正一下即可，call peak也很简单，后面的IGV可视化，DiffBind差异分析也能互相验证，再结合生物学的marker，对数据就基本掌控了。

基于这些solid的基础分析，就能做一些高级分析，主要是拿到peak、bam文件，就能在R里自由翱翔了。

 

目前M60这个项目的数据太丰富了，bulk的RNA-seq、ChIP-seq、ATAC-seq，Cut&Run，可重复性还行，缺点就是模型过于简单，来自几个cancer cell line。

 

注意：所有bulk的数据都要设定严格的read mapping filter，而single-cell则可以松一点，这取决于信噪比，bulk里存在大量的垃圾序列。

1. 去重复；
2. 去掉非unique的比对；
3. 去掉非paired的比对； 

1 去除低质量比对（MAPQ<30）
2 去除多重比对（一条read比对到基因组的多个位置）
3 去除PCR重复（不同reads比对到基因组的同一位置）
4 去除比对到黑名单区域的序列

RNA-seq一般不去重复（起始量高，扩增数少；某些RNA含量很高）
ChIP-seq一般去重复（起始量低）
call SNP一般去重复（PCR扩增错误影响SNP识别位点置信度）
PCR去重工具首选Picard
万事无绝对，还需参考起始量和PCR扩增数判断是否去重复。reads mapping覆盖均匀度可以判断是否需要去重复。
根源上解决去重复问题：起始量高，循环数少，reads能长不短，能双端不单端

 

ChIP-seq

本质就是研究DNA-protein互作，通过特异性抗体抓蛋白，对互作的DNA测序。

参考：

• hbctraining/Intro-to-ChIPseq
• ChIP-Seq上游数据处理一网打尽2023更新版

跟Cut&Run可以共享流程，主要差别就是rmDup、call peak时一定要用control，去掉背景噪音。

关于reads rm dup

照着教程里的分析，换成了samtools markdup，结果几乎所有read都被rm了，下次在花时间来测试这个功能。

换成picard之后一切都正常了，就是本来一行命令要换成两行。

 

ATAC-seq

本质是研究所有DNA开放区域，由Tn5酶切，线粒体完全裸露也会被捕获。

参考：

• 明码标价之ATAC-seq

这个目前还没正式跑过，scATAC-seq用的是cellranger的流程，用的是fragment来call peak，应该也是大同小异。

主要是多了一个BAM to tagAlign的步骤，然后还是用macs来call peak。

对于单端序列。直接用bed格式就可以；对于双端序列，推荐用bedpe格式。这两种格式都可以称之为tagAlign，可以作为macs的输入文件。 

 

关于bamtobed，需要sort by name！

*****WARNING: Query [...] is marked as paired, but its mate does not occur next to it in your BAM file. Skipping.

bedtools bamtobed -bedpe -i $sample.bam > $sample.bed

　　

 

参考：

• bam文件过滤和去重复