Calling Peaks on Integrated Multiome Dataset

[sylestiel]

@timoast
Have a couple of questions:

1. Should I be calling peaks on each dataset prior to integration or on the integrated dataset? This is with view to generating Coverage Plots and Peak to Gene Links.

2. When working with Mmul_10 genome (Rhesus Macaque) do I need to change the CallPeaks parameters in anyway?

Thanks,

[timoast]

You need to define a common set of peaks before combining your datasets. There are several ways you could do this, see here for an example.

For different species you can change the effective.genome.size parameter in CallPeaks()

[jamrute]

@timoast I had a question regarding integrating Multiome samples and peak calling as above. I followed the link you showed for merging files and I merged all our samples using the provided approach. Once we combine/QC as in that vignette is the idea then to following the following? If so, when creating the CreateChromatinAssay in the last step, we have to provide a path to fragment file, but since multiple samples were combined how do we go about it there?

Thanks

# call peaks using MACS2
peaks <- [CallPeaks](https://stuartlab.org/signac/reference/callpeaks)(pbmc)

# remove peaks on nonstandard chromosomes and in genomic blacklist regions
peaks <- [keepStandardChromosomes](https://rdrr.io/pkg/GenomeInfoDb/man/seqlevels-wrappers.html)(peaks, pruning.mode = "coarse")
peaks <- [subsetByOverlaps](https://rdrr.io/pkg/IRanges/man/findOverlaps-methods.html)(x = peaks, ranges = blacklist_hg38_unified, invert = TRUE)
# quantify counts in each peak
macs2_counts <- [FeatureMatrix](https://stuartlab.org/signac/reference/featurematrix)(
  fragments = [Fragments](https://stuartlab.org/signac/reference/fragments)(pbmc),
  features = peaks,
  cells = [colnames](https://rdrr.io/pkg/BiocGenerics/man/row_colnames.html)(pbmc)
)
# create a new assay using the MACS2 peak set and add it to the Seurat object
pbmc[["peaks"]] <- [CreateChromatinAssay](https://stuartlab.org/signac/reference/createchromatinassay)(
  counts = macs2_counts,
  fragments = fragpath,
  annotation = annotation
)