IDR分析

由于IDR分析比较严格，所以要是打算进行IDR分析的话，call peaks时候 设置阈值时候选择p值，且p < 0.01即可

#设置变量
sample_r1=SRR6322534
sample_r2=SRR6322535
sample=SRR6322534_SRR6322535
input=SRR6322538

#设置基因组大小
gsize=199000000

## call peaks for replicte 1
macs2 callpeak -t ./${sample_r1}.ff.bam -c ./$input.ff.bam -f BAM -n ${sample_r1} -g $gsize --keep-dup all -p 0.01

## call peaks for replicte 2
macs2 callpeak -t ./${sample_r2}.ff.bam -c ./$input.ff.bam -f BAM -n ${sample_r2} -g $gsize --keep-dup all -p 0.01

## call peaks for combined dataset，这里之间把两个bams 文件提供即可，不需要合并bams
macs2 callpeak -t ./${sample_r1}.ff.bam ./${sample_r2}.ff.bam -c ./$input.ff.bam -f BAM -n ${sample} -g $gsize --keep-dup all -p 0.01

## run idr
idr --samples ${sample_r1}_peaks.narrowPeak ${sample_r2}_peaks.narrowPeak --peak-list ${sample}_peaks.narrowPeak --input-file-type narrowPeak --output-file ./${sample}.idr --rank p.value --soft-idr-threshold 0.05 --plot

## get idr produced narrowPeak file，把结果转换成narrow peaks
cut -f 1-10 ${sample}.idr | sort -k1,1 -k2,2n -k3,3n >${sample}.idr.narrowPeak

## filter peaks against blacklist
bedtools intersect -a ${sample}.idr.narrowPeak -b ../ref/blacklist.bed -f 0.25 -v >tmp && mv tmp ${sample}.idr.narrowPeak

参考：基因课------表观基因组学之 ChIP-seq 数据分析