I am trying to compare two ChIP-seq datasets for histone modifications, however, the Spearman and Pearson correlations between two particular datasets does not correspond to their peak overlap (peak overlap is very low, while both pearson and spearman tests are showing strong correlation). I am performing this analysis on Usegalaxy, first using multibamSummary on both bam files followed by plotCorrelation. What could be the error here?
-
$\begingroup$ Can you please show any plots or data that have helped you come to the conclusion that there is a problem? $\endgroup$– gringer ♦Mar 1 at 18:49
-
1$\begingroup$ I can readily imagine a situation where the bulk data agree very well in terms of the overall genomic profile, but the statistical process by which you pick peaks yields different results. Knowing more about the samples would be handy- are they technical replicates, biological replicates, experimental arms? $\endgroup$– Maximilian PressMar 7 at 0:13
-
$\begingroup$ I agree with this, the correlation could be supported by the breadth of values across the data set. The actual r2-value is key: if it's 0.5 to 0.7 I wouldn't consider there to be a discrepancy. If it's >0.8 thats puzzling. At >0.9 there's definite selective sampling occurring because thats weird. The r alone (Spearman's) will be higher in most cases so I would use r2 as the guide. Summary a correlation can be stabilised by agreement between low and high values, so the exact correlation score is needed. $\endgroup$– M__ ♦Apr 19 at 15:58
1 Answer
Correlations like this between ChIP-seq datasets is common.
If you're calculating the correlation between BAM files across the entire genome, those are largely going to be the same (hence the need for input controls in ChIP-seq for performing the peak calling in the first place). Similarly, many ChIP-seq experiments only have 30 - 60% of reads within their peaks (FRiP), which means that 40 - 70% of reads don't fall within peak regions. This is going to heavily skew any correlation calculations that don't exclusively focus on the peaks themselves.
To more clearly compare ChIP-seq data against each other, global metrics like FRiP can help you diagnose measurement/sequencing issues. Looking at the correlation across ChIP-seq signal within peaks may be a more informative calculation than looking at all reads across the genome. Similarly, principal component analysis on signal within peaks can be helpful here.
DiffBind has some helper functions to perform these types of calculations that you can look at.