Many differentially accessible peaks identified, but not called as peak

[dmsalsgh97]

Hi, first, thanks for developing archR.
Its reproducible peak and TSS normalized bigwig strategy really work nicely for the case-control study I'm doing.

Now I'm analyzing a case-control study and I successfully clustered each cell type and divide it by conditions.
So I have Cell-Type_group divided clusters (ex: typeA_case and typeA_control)

I want to identify 'differentially accessible peaks' between the case-control groups, and I used getMarkerFeatures function with useGroups and bdgGroups options on PeakMatrix following the manual.

markerTest <- getMarkerFeatures(
    ArchRProj = ArchR_Snap, 
    useMatrix = "PeakMatrix",
    groupBy = "Snap_Celltype_by_condition",
    maxCells = 2000,
    binarize = TRUE,
    testMethod = "binomial",
    bias = c("TSSEnrichment", "log10(nFrags)"),
    useGroups = paste(cluster, "Dis", sep="_"), 
    bgdGroups = paste(cluster, "Normal", sep="_")
    )

For visualizing the output with IGV, I extracted my archR project's peak info (../PeakCalls/*-reproducilePeaks.gr.rds) to .bed files and coverage data to .bw files and load it on IGV.

But, When I searching differentially accessible peak regions on IGV, there exists "identified as differentially accessible peaks, but not called as peak"

As you can see on the upper IGV display, the lined DAR region is not peak called regions...

In this case, what is going wrong?
Thanks.

[rcorces]

You should review how the union reproducible peak set is genereated in ArchR
https://www.archrproject.com/bookdown/the-iterative-overlap-peak-merging-procedure.html

Because peak calls from multiple cell groups are merged, some peaks in the reproducible peak set will not exist in every cell type.

[dmsalsgh97]

Following your comment, I fully read the manual and paper again.
Can I extract DARs which peak called on the cluster using ArchR? ( for excluding that lined DAR case on the upper igv )

For example, if I have A_case cluster and A_control cluster,
I want to identify DARs between A_case and A_control that peak called on the A cluster ( A_case or A_control or both )
I want to do this because, in my think, DAR( potential causal regulatory elements ) should be called as peak on the cluster A.

Umm, In my think, using bedtools merge between A.bed and DAR.bed is a one possible option
Thanks!

[rcorces]

You can create whatever peakSet you want and add it to ArchR using addPeakSet()