Integration using Harmony as described in Drosophila eye-antennal disc preprint about cistopic

[cbravo93]

Hi @liaojinyue !

Very relevant question! There are 3 options:

1. Since I assume you have run the combined analysis (with the two conditions), you can try using the topic-cell matrix directly as input for Harmony (as if it was a PC matrix). We have applied this on several projects and are quite happy with the results.

2. Another option is to use the projection approach we used on the preprint, as you mention (see code below). Briefly, here we trained the model on 10X data, and then did a rough estimation of the topic contributions (trained on 10X) on the FACS cells (by multiplying the binary counts cell-region matrix by the region-topic contributions trained on 10X). Then we normalised the new contributions by Z-score, and combined to the Z-score normalised matrix topic-cell 10X matrix. Then this combined matrix you can use it as in 1); using as input for Harmony as if it was a PC matrix.

# Add cells to existing cisTopic model
## Make sure you use the same regions for both
EAD_10X <- readRDS('Your_path/cisTopicObject_EAD_10X.Rds')
EAD_Fluidigm <- readRDS('Your_path/cisTopicObject_EAD_Fluidigm.Rds')

## Select cells with at least 70% reads in peaks
FC1_binary <- [e-mail][,which([e-mail]$pct_ReadsInPeaks > 0.70)]

## Calculate region-topic distribution
topics <- [e-mail]$topics
beta <- 0.1 # Default beta
regiontopic_10X <- (topics + beta)/Matrix::rowSums(topics + beta)
colnames(regiontopic_10X) <- colnames(topics)
regiontopic_10X <- regiontopic_10X[,which(colnames(regiontopic_10X) %in% rownames(FC1_binary))]
FC1_binary <- F1C_binary[rownames(FC1_binary) %in% colnames(regiontopic_10X),]
FC1_binary <- F1C_binary[colnames(regiontopic_10X),]

## Calculate merged cell-topic distribution
celltopic_10X <- [e-mail]$document_sums
colnames(celltopic_10X) <- rownames([e-mail])
celltopic_FACS <- regiontopic_10X %*% F1C_binary
colnames(celltopic_FACS) <- colnames(F1C_binary)
rownames(celltopic_FACS) <- rownames(celltopic_10X)
combined_celltopic <- cbind(celltopic_10X, celltopic_FACS)
tech <- c(rep('10X', length(colnames(celltopic_10X))), rep('FACS', length(colnames(celltopic_FACS))))
names(tech) <- colnames(combined_celltopic)
tech <- as.data.frame(tech)

## Z-score and harmony
modelMat <- scale(combined_celltopic, center=TRUE, scale=TRUE)
combined_celltopic_harmony <- HarmonyMatrix(as.matrix(modelMat), tech, "tech", do_pca=FALSE)
colnames(combined_celltopic_harmony) <- colnames(combined_celltopic)
tSNE_coords <- Rtsne(t(as.matrix(combined_celltopic_harmony)), perplexity=100, seed=777, pca=FALSE)

PD: The rest of the code for the main analyses and to reproduce the main figures of the paper is available here. I will also add code for supplementary figures upon request.

3. Do you see any specific batch effect topic/s? Another option would be to remove this/these topic/s from the analysis.

Hope this is useful!

C

[liaojinyue]

Hi @cbravo93,

Thanks for your suggestion. I've tried option 1 and it works quite well.
Can you share with me why you chose option 2 in your preprint? Is there any advantage in doing integration this way?

Thanks,
Jason

[cbravo93]

Hi @liaojinyue

The main advantages of the option 2 is that (1) you don't need to rerun cisTopic on the combined data set and (2) you will work on the topic space of the reference data set, so if you have been working with those topics for a while makes it a bit easier. In the paper the goal was to project some FACS single-cell profiles in our existent model, more as a validation of the whole tissue model (to show that the FACS profiles cluster with the expected populations).

In this case, the main model already included cells matching the FACS cells, so the topics included information of these cell types. In your case, I guess you expect some change between your conditions, so if you project condition A into condition B your model won't include information about newly accessible regions and cell states on A (because they just weren't there in B).

To conclude, option 1 is really fast to check and can work if your reference model already includes cell types/states in your data sets to project (e.g. two runs in the same tissue); but if you expect changes between the reference and the others doing the combined analysis (followed by Harmony if you see batch effects) will give you more representative topics of what is going on.

Hope it is clear now!

C

[liaojinyue]

Hi @cbravo93,

Thank you for your detailed and clear explanation! That's extremely helpful. I will close this issue now.

Best,
Jason