ChipSeq: merge bam file before peak calling

Question

[Mike]

Hi all,

I have chipseq data for 6 samples (treatment/control), with replicate (so total 12 bam files). Now I am looking for peak calling using MACS2, I have doubt on handling replicates in macs2, so anyone please clear my doubt on handling replicates,

should I merge replicate before peak calling or analyse separately??

Thanks in advance.

Answer 1

[Devon Ryan]

I assume that your replicates are biological replicates. If they're technical then just merge the BAM files. Anyway, call peaks on each sample individually. You can then merge the peaks after filtering them and continue with whatever downstream processing is appropriate for your experiment.

[dariober]

Hi Devon,

You can then merge the peaks after filtering them

This is what I usually do but it opens up another can of worms about filtering peaks and defining a consensus set of peaks. One can call peaks at high stringency then get the union of all peaks, or lower stringency then accept peaks found in n out m files etc. (In practice I go for the solution the gives a consensus set that gives sensible results after some visual inspection with IGV).

[Mike]

Thanks Devon,

my replicates are technical, so I should just merge the BAM files.

Is this command Ok?

samtools merge  WTCHG_112.bam WTCHG_223245_112.bam WTCHG_223246_112.bam

[Ryan Dale]

FWIW, you can provide multiple BAMs to macs2 and it'll merge them for you, avoiding the extra step and extra storage.

[Mike]

Thanks Ryan, also looks good idea.

But how to use multiple files in macs2 command:

macs2 callpeak -t ChIP.bam -c Control.bam -f BAM -g hs -n test -B -q 0.01

[Ryan Dale]

Provide them as space-separated values. From macs2 callpeak -h:

Input files arguments:
  -t TFILE [TFILE ...], --treatment TFILE [TFILE ...]
                        ChIP-seq treatment file. If multiple files are given
                        as '-t A B C', then they will all be read and pooled
                        together. REQUIRED.
  -c [CFILE [CFILE ...]], --control [CFILE [CFILE ...]]
                        Control file. If multiple files are given as '-c A B
                        C', they will be pooled to estimate ChIP-seq
                        background noise.

[Devon Ryan]

Yup, that should work :)

Answer 2

[Constantine]

You can use a peakcaller that accepts biological replicates

I personally use PePr and I'm really happy with it. It has really high sensitivity and works very well both for sharp and broad regions

[Mike]

Thanks Constantine,

looks good, I will try..

[Mike]

Hi , I am using PePr , using following command but I got errors, could you please help me, what is wrong?

python PePr.py -c Chip_R1.bam -i Input_R1.bam –n R1_sample1 –f bam --peaktype=sharp

error:

Traceback (most recent call last):
  File "PePr.py", line 57, in <module>
    try: main(sys.argv)
  File "PePr.py", line 25, in main
    parameter, readData = optParser.process_opt(opt)
  File "optParser.py", line 84, in process_opt
    parameter = Parameters(opt)
  File "classDef.py", line 124, in __init__
    self.process_command_line_option(opt)
  File "classDef.py", line 142, in process_command_line_option
    self.file_format = opt.file_format.lower()
AttributeError: 'NoneType' object has no attribute 'lower'

[Constantine]

Sorry for the late reply

you need to parse replicates. For example

python PePr.py -c Chip_R1.bam,ChIP_R2.bam -i Input_R1.bam,Input_R2.bam –n R1_sample1 –f bam --peaktype=sharp