## ----echo=FALSE, results="hide", message=FALSE-------------------------------- knitr::opts_chunk$set(error=FALSE, message=FALSE, warning=FALSE) library(BiocStyle) ## ----------------------------------------------------------------------------- suppressMessages(library(scRepertoire)) ## ----eval=FALSE--------------------------------------------------------------- # S1 <- read.csv(".../Sample1/outs/filtered_contig_annotations.csv") # S2 <- read.csv(".../Sample2/outs/filtered_contig_annotations.csv") # S3 <- read.csv(".../Sample3/outs/filtered_contig_annotations.csv") # S4 <- read.csv(".../Sample4/outs/filtered_contig_annotations.csv") # # contig_list <- list(S1, S2, S3, S4) ## ----eval = F----------------------------------------------------------------- # contigs <- read.csv(".../outs/filtered_contig_annotations.csv") # # contig.list <- createHTOContigList(contigs, Seurat.Obj, group.by = "HTO_maxID") ## ----------------------------------------------------------------------------- data("contig_list") #the data built into scRepertoire head(contig_list[[1]]) ## ----eval=FALSE--------------------------------------------------------------- # for (i in seq_along(contig_list)) { # contig_list[[i]] <- stripBarcode(contig_list[[i]], column = 1, connector = "_", num_connects = 3) # } ## ----eval=FALSE--------------------------------------------------------------- # contig.output <- c("./Sample1/outs/", "./Sample2/outs/", "./Sample3/outs/") # contig.list <- loadContigs(dir = contig.output, # format = "TRUST4") ## ----------------------------------------------------------------------------- combined <- combineTCR(contig_list, samples = c("PY", "PY", "PX", "PX", "PZ","PZ"), ID = c("P", "T", "P", "T", "P", "T"), cells ="T-AB") ## ----------------------------------------------------------------------------- example <- addVariable(combined, name = "batch", variables = c("b1", "b1", "b2", "b2", "b2", "b2")) example[[1]][1:5,ncol(example[[1]])] # This is showing the first 5 values of the new column added ## ----------------------------------------------------------------------------- subset <- subsetContig(combined, name = "sample", variables = c("PX", "PY")) ## ----------------------------------------------------------------------------- quantContig(combined, cloneCall="strict", chain = "both", scale = TRUE) ## ----------------------------------------------------------------------------- quantContig_output <- quantContig(combined, cloneCall="strict", scale = TRUE, exportTable = TRUE) quantContig_output ## ----------------------------------------------------------------------------- quantContig(combined, cloneCall="gene", group.by = "ID", scale = TRUE) ## ----------------------------------------------------------------------------- abundanceContig(combined, cloneCall = "gene", scale = FALSE) abundanceContig(combined, cloneCall = "gene", group.by = "ID", scale = FALSE) ## ----------------------------------------------------------------------------- abundanceContig(combined, group.by = "ID", scale = TRUE) ## ----------------------------------------------------------------------------- lengthContig(combined, cloneCall="aa", chain = "both") lengthContig(combined, cloneCall="aa", chain = "TRA") ## ----------------------------------------------------------------------------- compareClonotypes(combined, numbers = 10, samples = c("PX_P", "PX_T"), cloneCall="aa", graph = "alluvial") ## ----------------------------------------------------------------------------- vizGenes(combined, gene = "V", chain = "TRB", plot = "bar", scale = TRUE) ## ----------------------------------------------------------------------------- #Peripheral Blood vizGenes(combined[c(1,3,5)], gene = "V", chain = "TRB", y.axis = "J", plot = "heatmap", scale = TRUE, order = "variance") #Tumor Infiltrating vizGenes(combined[c(2,4,6)], gene = "V", chain = "TRB", y.axis = "J", plot = "heatmap", scale = TRUE, order = "variance") ## ----------------------------------------------------------------------------- clonalHomeostasis(combined, cloneCall = "gene") clonalHomeostasis(combined, cloneCall = "aa") ## ----------------------------------------------------------------------------- clonalProportion(combined, cloneCall = "gene") clonalProportion(combined, cloneCall = "nt") ## ----------------------------------------------------------------------------- clonalOverlap(combined, cloneCall = "strict", method = "morisita") ## ----------------------------------------------------------------------------- clonesizeDistribution(combined, cloneCall = "strict", method="ward.D2") ## ----------------------------------------------------------------------------- clonalDiversity(combined, cloneCall = "gene", group.by = "sample", x.axis = "ID", n.boots = 100) ## ----------------------------------------------------------------------------- scatterClonotype(combined, cloneCall ="gene", x.axis = "PY_P", y.axis = "PY_T", dot.size = "total", graph = "proportion") ## ----------------------------------------------------------------------------- library(Seurat) library(scater) screp_example <- get(data("screp_example")) sce <- suppressMessages(UpdateSeuratObject(screp_example)) sce <- as.SingleCellExperiment(screp_example) #Seurat Format DimPlot(screp_example) ##Single Cell Experiment Format plotUMAP(sce, colour_by = "seurat_clusters") ## ----------------------------------------------------------------------------- table(screp_example$seurat_clusters) ## ----------------------------------------------------------------------------- screp_example <- combineExpression(combined, screp_example, cloneCall="gene", group.by = "sample", proportion = FALSE, cloneTypes=c(Single=1, Small=5, Medium=20, Large=100, Hyperexpanded=500)) sce <- combineExpression(combined, sce, cloneCall="gene", group.by = "sample", proportion = TRUE) ## ----eval=FALSE--------------------------------------------------------------- # #This is an example of the process, will not be evaluated during knit # TCR <- combineTCR(...) # BCR <- combineBCR(...) # list.receptors <- c(TCR, BCR) # # # seurat <- combineExpression(list.receptors, # seurat, # cloneCall="gene", # proportion = TRUE) ## ----------------------------------------------------------------------------- colorblind_vector <- colorRampPalette(rev(c("#0D0887FF", "#47039FFF", "#7301A8FF", "#9C179EFF", "#BD3786FF", "#D8576BFF", "#ED7953FF","#FA9E3BFF", "#FDC926FF", "#F0F921FF"))) DimPlot(screp_example, group.by = "Type") + NoLegend() + scale_color_manual(values=colorblind_vector(2)) ## ----------------------------------------------------------------------------- slot(screp_example, "meta.data")$cloneType <- factor(slot(screp_example, "meta.data")$cloneType, levels = c("Hyperexpanded (100 < X <= 500)", "Large (20 < X <= 100)", "Medium (5 < X <= 20)", "Small (1 < X <= 5)", "Single (0 < X <= 1)", NA)) DimPlot(screp_example, group.by = "cloneType") ## ----------------------------------------------------------------------------- clonalOverlay(screp_example, reduction = "umap", freq.cutpoint = 1, bins = 10, facet = "Patient") + guides(color = "none") ## ----------------------------------------------------------------------------- #ggraph needs to be loaded due to issues with ggplot library(ggraph) #No Identity filter clonalNetwork(screp_example, reduction = "umap", identity = "seurat_clusters", filter.clones = NULL, filter.identity = NULL, cloneCall = "aa") #Examining Cluster 2 only clonalNetwork(screp_example, reduction = "umap", identity = "ident", filter.identity = "C2", cloneCall = "aa") ## ----------------------------------------------------------------------------- screp_example <- highlightClonotypes(screp_example, cloneCall= "aa", sequence = c("CAVNGGSQGNLIF_CSAEREDTDTQYF", "NA_CATSATLRVVAEKLFF")) Seurat::DimPlot(screp_example, group.by = "highlight") ## ----------------------------------------------------------------------------- occupiedscRepertoire(screp_example, x.axis = "seurat_clusters") occupiedscRepertoire(screp_example, x.axis = "ident", proportion = TRUE, label = FALSE) ## ----------------------------------------------------------------------------- alluvialClonotypes(screp_example, cloneCall = "aa", y.axes = c("Patient", "ident", "Type"), color = "CAVNGGSQGNLIF_CSAEREDTDTQYF") + scale_fill_manual(values = c("grey", colorblind_vector(5)[5])) alluvialClonotypes(screp_example, cloneCall = "gene", y.axes = c("Patient", "ident", "Type"), color = "seurat_clusters") ## ----------------------------------------------------------------------------- library(circlize) library(scales) circles <- getCirclize(screp_example, group.by = "seurat_clusters") #Just assigning the normal colors to each cluster grid.cols <- hue_pal()(length(unique(screp_example$seurat_clusters))) names(grid.cols) <- unique(screp_example$seurat_clusters) #Graphing the chord diagram chordDiagram(circles, self.link = 1, grid.col = grid.cols) ## ----------------------------------------------------------------------------- subset <- subset(screp_example, Type == "T") circles <- getCirclize(subset, group.by = "ident") grid.cols <- scales::hue_pal()(length(unique(subset@active.ident))) names(grid.cols) <- levels(subset@active.ident) chordDiagram(circles, self.link = 1, grid.col = grid.cols, directional = 1, direction.type = "arrows", link.arr.type = "big.arrow") ## ----------------------------------------------------------------------------- sub_combined <- clusterTCR(combined[[2]], chain = "TRA", sequence = "aa", threshold = 0.85, group.by = NULL) ## ----------------------------------------------------------------------------- StartracDiversity(screp_example, type = "Type", sample = "Patient", by = "overall") ## ----------------------------------------------------------------------------- clonotypeBias(screp_example, cloneCall = "aa", split.by = "Type", group.by = "seurat_clusters", n.boots = 20, min.expand =1) ## ----------------------------------------------------------------------------- combined2 <- expression2List(screp_example, split.by = "seurat_clusters") combined3 <- expression2List(sce, split.by = "seurat_clusters") ## ----------------------------------------------------------------------------- clonalDiversity(combined2, cloneCall = "nt") clonalDiversity(combined3, cloneCall = "nt") ## ----------------------------------------------------------------------------- clonalHomeostasis(combined2, cloneCall = "nt") clonalHomeostasis(combined3, cloneCall = "nt") ## ----------------------------------------------------------------------------- clonalProportion(combined2, cloneCall = "nt") clonalProportion(combined3, cloneCall = "nt") ## ----------------------------------------------------------------------------- clonalOverlap(combined2, cloneCall="aa", method="overlap") clonalOverlap(combined3, cloneCall="aa", method="overlap") ## ----------------------------------------------------------------------------- sessionInfo()