If you want to find enriched peaks of drug treatment over the DMSO control then you would run a command like so:
macs2 bdgdiff --t1 TFdrug.bam --c1 TFdrug_input.bam --t2 TFDMSO.bam --c2 TFDMSO_input.bam --d1 --d2 -l <read lengths> --o-prefix TF
What you'll also want to do is specify --d1 & --d2 which is the reads filtered value from running macs2 predictd -i TFdrug.bam and macs2 predictd -i TFDMSO.bam and -l is your read lengths e.g. 36, 101, etc...
This macs2 bdgdiff should output 3 BEDs: one being t1 peaks enriched compared to t2, one being t2 peaks being enriched compared to t1, and one being peaks that are commonly enriched. You'll likely want to use the first BED file here.
Best of luck, Meer
• 0 views
•
link
written
by
meer.k.mustafa •
It may just be me, but I don't understand your post.
How can you use MACS2 to determine methylation? MACS2 is used to find enriched regions (peaks) in immunoprecipitated DNA, and does not check the methylation status of the bases, as far as I know.
To compare the immunoprecipitated regions in two ChIP-Seq experiments, without any background controls, I just run MACS2 callpeaks on both samples. I don't run MACS2 bdgdiff. So, just run MACS2 callpeak on both samples at the same time. You'll have to run it twice, inverting the treatment and the control, to get peaks enriched in both conditions.
But, this has nothing to do with methylation.
@ablanachetcohen My fault. I should mention that my libraries were generated with MeDIP-seq protocol.
Interesting. No, I guess it's just my ignorance. We always use bisulfite conversion, occasionally combined with ChIP-Seq (BisChIP-Seq). I hadn't heard of MeDIP.
My answer though would still be the same. Just run MACS2 callpeak twice, putting each sample as either treatment or control. If the background is the same in both samples, this approach is valid. I sometimes use this approach even when I have background controls for the samples, because MACS2 bdgdiff returns only a log10 probability ratio that is difficult to interpret. MACS2 callpeak returns much more detailed results that make it easier to subsequently filter the peaks.
I actually did think about using this approach, then I saw the tutorial of
bdgdiff. And I thought it may be the better tool to do the task. Anyway, thanks for the suggestion, I will try the "flipping peakcall" approach as well.