I'm trying to analyze a public dataset (GEO accession: GSE66318) which stems from a CHIP-Seq analysis by Vasudevan et al. 2015. Specifically, I'm interested in the peak-calling results they provide for the acetylation marks for both treated and control. These results are in the form of BED files which contain coordinates for high-confidence peaks for each sample and each replicate, and two bigWig files per replicate containing the (I suppose high-confidence) peaks called by MACS2. The two bigWig files per replicate stem from the fact that apparently MACS2 doesn't do background subtraction nor normalization on the results, so we have a background peaks file (lambda) and a CHIP-Seq peaks file ("treatment") for each replicate.
What I want to do is analyze this dataset to try and find out if there are significant differences in the acetylation levels between the NO treated and untreated groups for a particular gene, after background-subtracting and normalizing each replicate and taking in account the biological replicates for each group. Is it possible to do that with only the bigWig and BED files in hand? If so, can you give me some direction? Some of the steps for account for biological replicates that I've seen by following this tutorial assume I have the original bam alignment files. The same for background-filtering the peak-calling results. I've been trying to explore the deepTools toolset to see if I can generate the results I want from the files I have, but I'm feeling a bit lost.
Thanks in advance.
