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)
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")