library(Seurat) library(tibble) library(magrittr) # Read in data data <- data.table::fread("~/Downloads/~.txt.gz", data.table = F) # Gene names are in the first column so we need to move them to rownames data <- data %>% column_to_rownames("GENE") # Create Seurat NT1<-CreateSeuratObject(counts = data, min.cells = 3, min.features = 350, project = "NT1") data2 <- data.table::fread("C:/Users/~/Downloads/~.txt.gz", data.table = F) data2 <- data2 %>% column_to_rownames("GENE") NT2<-CreateSeuratObject(counts = data2, min.cells = 3, min.features = 350, project = "PS") data3 <- data.table::fread("C:/Users/~/Downloads/~.txt.gz", data.table = F) data3 <- data3 %>% column_to_rownames("GENE") TT1<-CreateSeuratObject(counts = data3, min.cells = 3, min.features = 350, project = "PS") data4 <- data.table::fread("C:/Users/~/Downloads/~.txt.gz", data.table = F) data4 <- data4 %>% column_to_rownames("GENE") TT2<-CreateSeuratObject(counts = data4, min.cells = 3, min.features = 350, project = "PS") NT1[["percent.mt"]] <- PercentageFeatureSet(NT1, pattern = "^MT-") NT2[["percent.mt"]] <- PercentageFeatureSet(NT2, pattern = "^MT-") TT1[["percent.mt"]] <- PercentageFeatureSet(TT1, pattern = "^MT-") TT2[["percent.mt"]] <- PercentageFeatureSet(TT2, pattern = "^MT-") VlnPlot(NT1, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3) VlnPlot(NT2, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3) VlnPlot(TT1, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3) VlnPlot(TT2, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3) NT1 <- subset(NT1, subset = nFeature_RNA > 100 & nFeature_RNA < 4000 & percent.mt < 20) NT2 <- subset(NT2, subset = nFeature_RNA > 100 & nFeature_RNA < 4000 & percent.mt < 20) TT1 <- subset(TT1, subset = nFeature_RNA > 100 & nFeature_RNA < 4000 & percent.mt < 20) TT2 <- subset(TT2, subset = nFeature_RNA > 100 & nFeature_RNA < 4000 & percent.mt < 20) NT1@meta.data$condition <- c('NT1') NT2@meta.data$condition <- c('NT2') TT1@meta.data$condition <- c('TT1') TT2@meta.data$condition <- c('TT2') dat <- merge(NT1, y = c(NT2, TT1, TT2), add.cell.ids = c("NT1", "NT2", "TT1", "TT2"), project = "PS") PS.list<-SplitObject(dat, split.by = "condition") PS.list <- lapply(X = PS.list, FUN = function(x) { x <- NormalizeData(x) x <- FindVariableFeatures(x, selection.method = "vst", nfeatures = 4000) }) features <- SelectIntegrationFeatures(object.list = PS.list) PS.anchors <- FindIntegrationAnchors(object.list = PS.list, anchor.features = features) PS.combined <- IntegrateData(anchorset = PS.anchors) DefaultAssay(PS.combined) <- "integrated" # Run the standard workflow for visualization and clustering PS.combined<- ScaleData(PS.combined, verbose = FALSE) PS.combined<- RunPCA(PS.combined, npcs = 30, verbose = FALSE) PS.combined<- RunUMAP(PS.combined, reduction = "pca", dims = 1:30) PS.combined<- FindNeighbors(PS.combined, reduction = "pca", dims = 1:30) PS.combined<- FindClusters(PS.combined, resolution = 0.5) p1 <- DimPlot(PS.combined, reduction = "umap", group.by = "condition") p2 <- DimPlot(PS.combined, reduction = "umap", label = TRUE, repel = TRUE) p1 + p2 DimPlot(PS.combined, reduction = "umap", split.by = "condition") DefaultAssay(PS.combined) <- "RNA" seurat <- NormalizeData(object = PS.combined, normalization.method = "LogNormalize", scale.factor = 1e4) seurat <- FindVariableFeatures(seurat, selection.method = "vst", nfeatures = 4000) top10 <- head(VariableFeatures(seurat), 10) [1] "MYH11" "LTF" "OLFM4" "SPARCL1" "APOE" "SPP1" "DCN" "GPNMB" "RGS5" "ACTA2" [1] "MMP7" "OLFM4" "MYH11" "CPA3" "LTF" "IL1B" "SPARCL1" "GPNMB" "KIT" "SCGB1A1" nk.markers <- FindConservedMarkers(PS.combined, ident.1 = 6, grouping.var = "condition", verbose = FALSE) head(nk.markers) library(ggplot2) library(cowplot) theme_set(theme_cowplot()) plots <- VlnPlot(PS.combined, features = c("KRT8", "KRT18", "CD24"), split.by = "condition", group.by = "celltype", pt.size = 0, combine = FALSE) wrap_plots(plots = plots, ncol = 1) plots <- VlnPlot(PS.combined, features = c("ACTA2", "COL1A1", "KRT19", "PIK3C2A"), split.by = "condition", idents = "9", pt.size = 0, combine = FALSE) FeaturePlot(BC.combined, features = c("COL1A1", "COL1A2", "ACTA2", "THY1", "EPCAM", "KRT8", "KRT18", "CD24", "PTPRC"), min.cutoff = "q9") > FeaturePlot(BC.combined, features = c("COMMD1", "NFKB1", "HIF1A", "SLC16A1", "SLC16A3", "TIGAR", "ZEB1", "CDH1", "CAV1"), min.cutoff = "q9") markers.to.plot <- c("COL1A1", "COL1A2", "PDPN", "ACTA2", "EPCAM", "CDH1", "MSMB", "CAV1", "FTH1", "RPLP1", "HSPB1", "KLK3", "KRT5", "LCN2", "PSCA", "PTPRC", "CD68", "CD163", "CD14", "SLC16A3", "ITGAX", "CD3D", "CD8A", "MS4A1", "CD79A", "VWF", "CD200") DotPlot(PS.combined, features = markers.to.plot, cols = c("blue", "red", "green", "#00BFC4"), dot.scale = 8, split.by = "condition") + RotatedAxis() markers.to.plot <- c("COL1A1", "COL1A2", "PDPN", "ACTA2", "EPCAM", "CDH1", "MSMB", "CAV1", "FTH1", "RPLP1", "HSPB1", "KLK3", "KRT19", "PKM", "OS9", "KLF4", "KLF9", "PIK3C2A", "PTPRC", "CD68", "VWF", "CD200") > DotPlot(PS.combined, features = markers.to.plot, cols = c("blue", "red", "green", "#00BFC4"), dot.scale = 8, split.by = "condition") + RotatedAxis() C4 <- read.csv("C:/Users/~/Desktop/gsm5353232TC/nk.markers4.csv") CC4 <- merge(NTF, C4) write.csv(CC4, "cluster4.csv") |