CHIP-Seq(6)：查看两次重复中的peak情况

1.去除.narrowPeak文件首行

sed -i '1d' WT_rep1_peaks.narrowPeak 
sed -i '1d' WT_rep2_peaks.narrowPeak 

2.使用bedtools粗略估计两次重复的peak交集,重叠1bp就算有交集

bedtools intersect -a WT_rep1_peaks.narrowPeak -b WT_rep2_peaks.narrowPeak -wo > overlaps.bed

3.粗略估计后可以画韦恩图进行可视化

intervene venn -i *.narrowPeak --output venn --save-overlaps

image.png

4.对MACS3的结果文件narrowPeak根据-log10(p-value)进行排序

sort -k8,8nr WT_rep1_peaks.narrowPeak >WT_rep1_peaks.narrowPeak
sort -k8,8nr WT_rep2_peaks.narrowPeak >WT_rep2_peaks.narrowPeak

5.对生物学重复样本间的peak进行鉴定,查看两次重复的peak的IDR(不可重复率)

此次实验是具有生物学重复样本，处理前需要对重复样本的共有peak进行鉴定，用IDR的方法获得高重复性的peaks

idr --samples WT_rep1_peaks.narrowPeak WT_rep2_peaks.narrowPeak  --output-file idr_peak.narrowPeak --rank p.value --plot  --idr-threshold 0.05 --log-output-file sample.idr.log
#如果想看IDR<0.05的，可以通过第5列信息过滤：
awk '{if($5 >= 540) print $0}' idr_peak.narrowPeak | wc -l 

输出文件解读：
输出文件包括：

idr_peak.narrowPeak
sample-idr.log
sample-idr.png
（1）idr_peak.narrowPeak
idr_peak.narrowPeak是common peaks的结果输出文件，格式与输入文件格式类似，只是多了几列信息。前10列是标准的narrowPeak格式文件，包含重复样本整合后的peaks信息。
第5列：包含缩放的 IDR 值
score int
如min(int(log2(-125IDR), 1000)，那么IDR=0，缩放的IDR就是1000；IDR=0.05, int(-125log2(0.05)) = 540；IDR=1.0, int(-125log2(1.0) = 0。
<meta charset="utf-8">

• 11列和第12列：分别是local和global IDR值
  • col11: localIDR float -log10(Local IDR value)
  • col12: globalIDR float -log10(Global IDR value)
    global IDR值是第5列中用于计算缩放的IDR值，类似于对p值进行多个假设校正以计算FDR；local IDR类似于属于不可重复噪声部分的峰值的后验概率。
• 第13、14、15、16列：是重复样本1相关的信息。
  • col13: rep1_chromStart int
  • col14: rep1_chromEnd int
  • col15: rep1_signalValue float
  • col16: rep1_summit int
• 第17-20列：是重复样本2相关的信息，和13-16四列信息类似。

其他列信息如下：

• col1: chrom string
• col2: chromStart int
• col3: chromEnd int
• col4: name string
• col6: strand [+-.] Use '.' if no strand is assigned.
• col7: signalValue float
• col8: p-value float
• col9: q-value float
• col10: summit int

wc -l *-idr 计算下common peaks的个数，接着可再计算下与总peaks的比率。

2）sample-idr.log
log文件会给出peaks通过IDR < 0.05的比率

cat sample.idr.log 
Initial parameter values: [0.10 1.00 0.20 0.50]
Final parameter values: [1.22 1.01 0.89 0.29]
Number of reported peaks - 37/299 (12.4%)

Number of peaks passing IDR cutoff of 0.05 - 37/299 (12.4%)

(3）sample-idr.png
png文件包括4个图
左上： Rep1 peak ranks vs Rep2 peak ranks, 没有通过特定IDR阈值的peaks显示为红色。
右上：Rep1 log10 peak scores vs Rep2 log10 peak scores，没有通过特定IDR阈值的peaks显示为红色。
下面两个图： Peak rank vs IDR scores，箱线图展示了IDR值的分布，默认情况下，IDR值的阈值为-1E-6。

image.png