Reproducing marker gene visualizations

[koenvandenberge]

Dear Developers,

We have been using ArchR to interpret our scATAC-seq dataset. The marker genes identified by ArchR using the gene activity scores make a lot of sense; however, I am unable to reproduce the visualizations of the marker genes as provided by ArchR when trying to implement them myself.

The default, as recommended by your documentation, is to do

p <- plotEmbedding(
    ArchRProj = archProj, 
    colorBy = "GeneScoreMatrix", 
    name = markerGenes, 
    embedding = "UMAP_cor",
    quantCut = c(0.01, 0.95),
    size=1/3,
    plotAs="points"
)
# clean
p2 <- lapply(p, function(x){
    x + guides(color = FALSE, fill = FALSE) + 
    theme_ArchR(baseSize = 6.5) +
    theme(plot.margin = unit(c(0, 0, 0, 0), "cm")) +
    theme(
        axis.text.x=element_blank(), 
        axis.ticks.x=element_blank(), 
        axis.text.y=element_blank(), 
        axis.ticks.y=element_blank()
    )
})
cowplot::plot_grid(plotlist=p2)

These indeed look nice and the markers make sense.
I then tried to reproduce the visualization such that I would have more flexibility in tweaking them. I used standard

## make plots manually
umapEmbedding <- getEmbedding(archProj, embedding="UMAP_cor")
colnames(umapEmbedding) <- paste0("UMAP",1:2)
featureDF <- ArchR:::.getFeatureDF(head(ArrowFiles, 2), "GeneScoreMatrix")
featureDF <- featureDF[which(featureDF$name %in% markerGenes),]
geneActMat <- ArchR:::.getPartialMatrix(ArrowFiles, 
            featureDF = featureDF, threads = 1, useMatrix = "GeneScoreMatrix", 
            cellNames = rownames(archProj@cellColData), progress = FALSE)


  curDF <- umapEmbedding
  curDF$markerGene <- log1p(geneActMat[1,])
  curP <- ggplot(curDF, aes(x=UMAP1,
                            y=UMAP2,
                            col=markerGene)) +
    theme_classic() +
    ggtitle(markerGenes[1]) +
    geom_point(pch=16) +
    scale_color_gradient( breaks=c(0, quantile(curDF$markerGene, probs=seq(0.05, 0.95, length=98)) , 
max(curDF$markerGene))) +
    theme(legend.position = "none")

This, however, is giving me a radically different picture, while I thought that the visualization should be the same. What am I missing in order to be able to reproduce these figures?