ChIP-seq分析的一般流程方法
现在ChIP-seq的数据基本是最常见的测序数据类型之一，主要有Transcription factor ChIP-seq和Histone ChIP-seq。前者是看转录因子的结合位置，后者是组蛋白修饰发生的位置。下面分享一下一般流程。
ChIP-seq
- 质控 (quality control)
首先要看一下ChIP-seq数据的质量，数据的信号最好比background要强很很多。一般要有control，这样call peaks更准确可信， control主要有Input DNA 和 IgG两种，前一种更常用。
检测质量的一些方式：
1). peaks中reads的数量，如果peaks的reads普遍较少，则质量一般。
2). peaks信号高，背景低。
3). 测序深度深 。
4). Diverse library (与重复duplications有关，如下图)
library
4). 有重复并且与重复之间相似性较高…
……
做质控的软件方法：
1). ChIPQC (T Carroll, Front Genet, 2014.)
2). SPP package - Unix/Linux (PV Karchenko, Nature Biotechnol, 2008.)
3). ENCODE中的标准流程
- 序列比对 (mapping of fastq)
序列比对一般用BWA或者Bowtie2，两者效果差不多。BWA的bwa samse（单端数据）和bwa sampe (双端数据) 跑的速度比较慢，但是效果很不错，用法如下：
bwa index reference.fa # 建立索引 -p可设置前缀，不设置前缀就是reference.fa。
# 单端数据
bwa aln -t 8 reference.fa test.fq.gz > test.sai
bwa samse -n 10 reference.fa test.sai test.fq.gz > test_se.sam
# 双端数据：
bwa aln reference.fa test_reads1.fq > test1.sai
bwa aln reference.fa test_reads2.fq > test2.sai
bwa sampe reference.fa test1.sai test2.sai test_reads1.fq test_reads2.fq > test_pe.sam
BWA的mem，速度很快：
bwa mem reference.fa reads.fq > test_se.sam # 单端
bwa mem reference.fa read1.fq read2.fq > test_pe.sam # 双端
bowtie2的用法：
bowtie2-build reference.fa index # 创建index
bowtie2 -p 8 -x index -1 test_read1.fq -2 test_read2.fq -S test_se.sam # 单端比对
bowtie2 -p 8 -x index -1 test_read1.fq -2 test_read2.fq -S test_pe.sam # 双端比对
效果个人感觉差不多。
- 去除重复 (remove duplicates)
由于PCR实验存在不可避免的实验误差，所以会存在重复 (duplicates)。比如两条不同的reads，起止位置完全一致。比如：
example_dup
其中第二条已经被picard标注出来了。被标注的第二列flag会加1024。
去重的软件中samtools rmdup （基本已不用），samtools markdup（更新后的）和picard最常用。rmdup效果不怎么好，而且如果有遇到相同位置的reads， 会优先选择质量高的那一条read。picard与samtools markdup效果相似（仿佛调用的同一个？并不确定）。都可以标记重复，也可以选择直接去掉。以下是用法：
samtools markdup -@ 8 -r test.bam filter_test.bam # -r是直接去掉重复，不加是直接标记
picard去重有三种方式可选，在DUPLICATE_SCORING_STRATEGY参数中，分别是SUM_OF_BASE_QUALITIES, TOTAL_MAPPED_REFERENCE_LENGTH和RANDOM。即当有重复时分别选择留下总碱基质量最高的、匹配上参考基因组最长的和随机。
picard MarkDuplicates I=test.bam O= filter_test.bam M=dup_metrics.txt REMOVE_DUPLICATES=true
- peak calling
peaks是reads信号比较强的区域，也就是我们找到的转录因子或者组蛋白修饰最有可能结合的地方。call peaks仍然有不少软件，比较常用的是MACS2和Hotspot2。
示例：
macs2 callpeak -t test.bam -c control.bam -f BAM -g hs -n test -B -q 0.01
针对不同的数据考虑用不同的参数。
- 下游分析 (downstream analysis)
分析完之后下游可以做的事情很多，视情况而定。可以同时分析DNase-seq或者ATAC-seq的数据，看转录因子与染色质开放区的关系；或者Homer等工具注释peaks，看不同转录因子/组蛋白修饰之间的关系，或者分析TF的target gene。也可以用MEME等做motif分析。
欢迎关注公众号！
生信编程日常
