== Demultiplexing hashtag oligos (HTOs) == === Demultiplexing with hashtag oligos (HTOs) and analysis with Seurat === Follow the the Seurat v3 vignette on hashing[[https://satijalab.org/seurat/archive/v3.1/hashing_vignette.html]], but on the demultiplex step use kfunc = “kmeans”. These are the main steps: {{{ library(Seurat) }}} Read the data {{{ L13_892.data <- Read10X(data.dir = "/PathTo/L13_892_HTO/outs/filtered_feature_bc_matrix") }}} Explore the data {{{ head(L13_892.data$`Antibody Capture`@Dimnames[1]) rownames(x = L13_892.data[["Antibody Capture"]]) }}} Create Seurat object with the gene expression data {{{ L13_892 <- CreateSeuratObject(counts = L13_892.data[["Gene Expression"]]) }}} Normalize RNA data with log normalization {{{ L13_892 <- NormalizeData(L13_892) }}} Add the Percentage of mitochondrial reads {{{ L13_892[["percent.mt"]] <- PercentageFeatureSet(object = L13_892, pattern = "^mt-") }}} Add HTO data as a new assay independent from RNA {{{ L13_892 <- NormalizeData(L13_892, assay = "HTO", normalization.method = "CLR") }}} Demultiplex {{{ L13_892 <- HTODemux(L13_892, assay = "HTO", positive.quantile = 0.99, kfunc = "kmeans") }}} Make ridge plots to see how well the tags are separating {{{ Indents_Demux_892 <- Idents(L13_892) pdf("./RindgePlots_892.pdf", w=11, h=8.5) RidgePlot(L13_892, assay = "HTO", features = rownames(L13_892[["HTO"]])[1:2],ncol = 2) RidgePlot(L13_892, assay = "HTO", features = rownames(L13_892[["HTO"]])[3:4],ncol = 2) # Group cells based on the max HTO signal #This is just a quality control. # This is not the way of separating cells by barcode Idents(L13_892) <- "HTO_maxID" RidgePlot(L13_892, assay = "HTO", features = rownames(L13_892[["HTO"]])[1:2],ncol = 2) RidgePlot(L13_892, assay = "HTO", features = rownames(L13_892[["HTO"]])[3:4],ncol = 2) dev.off() Idents(L13_892) <- Indents_Demux_892 }}} See how many cells of each type do I have {{{ table(Idents(object = L13_892)) }}} Now that we have separated the different cell populations we would continue with the standard Seurat Analysis. If we are combining several lanes and want to keep track of which negatives and doublets are coming from each lane, we may want to rename the cell IDs before merging several objects: {{{ table(Idents(object = L13_892)) new.cluster.ids <- c("Doublet892", "BFUEs-24h-PlusDex", "BFUEs-24h-PlusGalun", "BFUEs-24h-PlusTGFb", "BFUEs-24h-RegMed", "Negative892") names(new.cluster.ids) <- levels(L13_892) L13_892 <- RenameIdents(L13_892, new.cluster.ids) table(Idents(object = L13_892)) head(Idents(L13_894)) table(Idents(object = L13_894)) new.cluster.ids <- c("Doublet894", "Condition1", "Condition2", "Condition3", "Condition4", "Negative894") names(new.cluster.ids) <- levels(L13_894) L13_894 <- RenameIdents(L13_894, new.cluster.ids) table(Idents(object = L13_894)) }}} Merge all the sequencing lanes (the other samples 894, 896, 898 and 900 have been processed the same way than 892) {{{ all <- merge(x = L13_892, y = c(L13_894,L13_896,L13_898,L13_900)) }}}