Question: How is FDR calculated in MACS?
4
FreeMindAlex100 wrote:

Hello All,

My question concerns the calculation of empirical FDR by MACS peak caller. According to the authors "For a particular P value threshold, the empirical FDR is then calculated as the number of control peaks passing the threshold divided by the number of ChIP-seq peaks passing the same threshold". I've got a problem grasping this. How can there be multiple peaks for a given P-value (estimated as significant enrichment above background) if P-value is supposed to refer to a single peak whether in control or ChIP sample?

Thanks

chip-seq • 5.8k views
modified 3.2 years ago by igor8.8k • written 4.4 years ago by FreeMindAlex100
10
Gary470 wrote:

Hi,

I have a similar question before. Actually, the sentences you show are not easy to understand for a new hand. The below are better descriptions from their papers and an example I try to explain for your reference.

Feng, J. X., T. Liu, B. Qin, Y. Zhang & X. S. Liu (2012) Identifying ChIP-seq enrichment using MACS. Nat. Protoc., 7, 1728.

Estimating the empirical false discovery rate by exchanging ChIP-seq and control samples. When a control sample is available, MACS can also estimate an empirical FDR for every peak by exchanging the ChIP-seq and control samples and identifying peaks in the control sample using the same set of parameters used for the ChIP-seq sample. Because the control sample should not exhibit read enrichment, any such peaks found by MACS can be regarded as false positives. For a particular P value threshold, the empirical FDR is then calculated as the number of control peaks passing the threshold divided by the number of ChIP-seq peaks passing the same threshold.

Zhang, Y., T. Liu, C. A. Meyer, J. Eeckhoute, D. S. Johnson, B. E. Bernstein, C. Nussbaum, R. M. Myers, M. Brown, W. Li & X. S. Liu (2008) Model-based Analysis of ChIP-Seq (MACS). Genome Biol., 9.

For a ChIP-Seq experiment with controls, MACS empirically estimates the false discovery rate (FDR) for each detected peak using the same procedure employed in the previous ChIP-chip peak finders MAT  and MA2C . At each p- value, MACS uses the same parameters to find ChIP peaks over control and control peaks over ChIP (that is, a sample swap). The empirical FDR is defined as Number of control peaks / Number of ChIP peaks…. To compare their prediction specificity, we swapped the ChIP and control samples, and calculated the FDR of each algorithm as Number of control peaks / Number of ChIP peaks using the same parameters for ChIP and control.

An example

For example, MACS uses a H3K27ac ChIP sample as the treat and an Input sample as the control to identify a peak, called as peak X, with p-value = 0.00024 based on the parameters you set. Overall, there are 1,000 peaks whose p-value ≦ 0.00024 (i.e.multiple peaks for a given P-value). After that, MACS uses the Input sample as the treat and the H3K27ac ChIP sample as the control to identify peaks again. If totally there are 48 peaks in the Input over the H3K27ac ChIP whose p-value ≦ 0.00024. The FDR for the peak X = 48 / 1000 = 0.048. It means that based on the same threshold, only 4.8 peaks out of 100 peaks in our H3K27ac ChIP sample are false positive, if we believe that all Input peaks are not real enrichment.

Thanks Gary

1
igor8.8k wrote:

I would add that FDR calculation is different depending on MACS version. According to the author:

In MACSv1.4, FDR is calculated by swapping treatment and control. MACSv14 assumes all peaks called under a cutoff in this way are false positives. So it can calculate empirical FDR for each pvalue cutoff. However this method is hugely influenced by unbalanced sequencing depth. For example, if control sample is much larger than treatment, FDR would be overestimated so you would have less 'good' peaks above cutoff.

In MACSv2, I use another approach. First, pvalue are calculated at every basepair in the genome, then I adopt Benjamini-Hochberg to correct multiple comparisons, convert pvalue into qvalue or minimum FDR that the peak is significant. This method is more robust in my experience.