In ChIP-Seq Drosophila data, there is a large difference between the Input samples. Calling peaks when one is taken as the treatment, and another as the control yields thousands of peaks.
Plotting the bam files in Ngsplot reveals that after normalization, in some of the Inputs there is much more binding in the genebody area than in others.
The inputs are not from exactly the same condition, but from very similar ones (if possible, I'd rather not go into details). In theory, there should be no large differences between the conditions.
How can I determine if this is a biological difference, or some technical artifact / problem with the analysis?
The data is 43bp paired-end; it was aligned in bowtie -m mode (yields uniquely mapped reads) to BDGP6.28. Bam files plotted with Ngsplot with dm6 genome (it's an older release than BDGP6.28 -Ensembl 79). Fragment size in ngsplot was set to 250, which was roughly the fragment size calculated by macs2. Blacklisted areas were not removed from the bam file. In this specific image duplicate reads were removed, but removing them makes no difference.
(Another issue in the plots is that the binding demonstrates peaks in the TSS, and the TES. As far as I understand, open chromatin is more easily sonicated than closed chromatin, and hence the binding is not expected to be even across the genebody. Here is an article showing that)