Thank you for your inputs. These variants do not always co-occupy a region, rather, at TSS of genes that are going to be differentially expressed across developmental timepoints. We can safely say that they preferentially bind transcriptionally active regions and make the nucleosome unstable.
Coming to "capturing co-ocupancy", it is basically inspired by Fig 1e in this paper: https://doi.org/10.1038/s41467-025-57719-4 where they showed bivalent regions using a heatmap. What we have done is basically take a union of peakfiles of both variants and then we computed a matrix of their signal enrichment (normalized bigwig tracks) on the peakfile. The core idea was, when plotted on an enrichment heatmap, the peaks that are overlapping (have co-occupancy) should form a seperate kmean cluster while the exclusive ones will form other. The results were not as expected because one variant has a much broader distribution and signal intensity than the other.
I have not heard about the method you described, can you please share link to any study that used this method?