I have performed peak calling on a number of separate ChIP-seq experiments and would like to harmonise the peaks from each of the experiments in order to convert my data into a matrix for further analysis. I would like to do this by taking the union of the called peaks (available in BED / narrowPeak / MACS2 format) and retrieving the enrichment score / q-value for all of the experiments on each peak. Is this possible without writing my own extension to MACS2?
$\begingroup$
$\endgroup$
3
asked Nov 26, 2018 at 5:13
-
1$\begingroup$ This is going to be an custom script so you're not going to touch MACS2 itself. The only question is how you want to summarize q-values from different samples. Other than that this is a job for bedtools/bedops. $\endgroup$– Devon RyanCommented Nov 26, 2018 at 11:34
-
$\begingroup$ Here's why I think I still need MACS - say Sample A has a peak at chrN:100-200, but Sample B does not. In order to avoid missing values, I would like to get the (not significant) enrichment score from that peak for Sample B using the MACS2 algorithm. $\endgroup$– Scott GiganteCommented Nov 27, 2018 at 14:51
-
$\begingroup$ Ah, if you really want that then yes you'll need to use the (undocumented) MACS2 API. $\endgroup$– Devon RyanCommented Nov 27, 2018 at 15:04
Add a comment
|
1 Answer
$\begingroup$
$\endgroup$
You could try setting a high p-value threshold when you're calling peaks to retain the "non-significant" values. Something like:
macs2 callpeak -t <ChIP>.bam -c <Control>.bam -f BED -g hs -n test -B -p 0.7
You probably don't want to set the threshold to 1, since you'll likely want to avoid some noise.
But you can then use these "non-significant" values in combination with a bedtools map command, to map them to your predefined union of "actual" peaks.
