comparison Seurat.R @ 2:321bdd834266 draft

"planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/seurat commit 3cf715ec11e2c9944f46572e324e5b2db5aa151f"
author iuc
date Thu, 19 Dec 2019 02:42:56 -0500
parents 7319f83ae734
children 764f076e9d52
comparison
equal deleted inserted replaced
1:7319f83ae734 2:321bdd834266
21 #' tsne: "" 21 #' tsne: ""
22 #' heatmaps: "" 22 #' heatmaps: ""
23 #' --- 23 #' ---
24 24
25 #+ echo=F, warning = F, message=F 25 #+ echo=F, warning = F, message=F
26 options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } ) 26 options(show.error.messages = F, error = function(){cat(geterrmessage(), file = stderr()); q("no", 1, F)})
27 showcode <- as.logical(params$showcode) 27 showcode <- as.logical(params$showcode)
28 warn <- as.logical(params$warn) 28 warn <- as.logical(params$warn)
29 varstate <- as.logical(params$varstate) 29 varstate <- as.logical(params$varstate)
30 vlnfeat <- as.logical(params$vlnfeat) 30 vlnfeat <- as.logical(params$vlnfeat)
31 featplot <- as.logical(params$featplot) 31 featplot <- as.logical(params$featplot)
32 PCplots <- as.logical(params$PCplots) 32 PCplots <- as.logical(params$PCplots)
33 tsne <- as.logical(params$tsne) 33 tsne <- as.logical(params$tsne)
34 heatmaps <- as.logical(params$heatmaps) 34 heatmaps <- as.logical(params$heatmaps)
35 35
36 # we need that to not crash galaxy with an UTF8 error on German LC settings. 36 # we need that to not crash Galaxy with an UTF-8 error on German LC settings.
37 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") 37 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
38 38
39 39
40 #+ echo = F, warning = `warn`, include =`varstate` 40 #+ echo = F, warning = `warn`, include =`varstate`
41 min_cells <- as.integer(params$min_cells) 41 min_cells <- as.integer(params$min_cells)
55 print(paste0("Resolution: ", resolution)) 55 print(paste0("Resolution: ", resolution))
56 print(paste0("Minimum percent of cells", min_pct)) 56 print(paste0("Minimum percent of cells", min_pct))
57 print(paste0("Logfold change threshold", logfc_threshold)) 57 print(paste0("Logfold change threshold", logfc_threshold))
58 58
59 #+ echo = FALSE 59 #+ echo = FALSE
60 if(showcode == TRUE){print("Read in data, generate inital Seurat object")} 60 if(showcode == TRUE) print("Read in data, generate inital Seurat object")
61 #+ echo = `showcode`, warning = `warn`, message = F 61 #+ echo = `showcode`, warning = `warn`, message = F
62 counts <- read.delim(params$counts, row.names=1) 62 counts <- read.delim(params$counts, row.names = 1)
63 seuset <- Seurat::CreateSeuratObject(counts = counts, min.cells = min_cells, min.features = min_genes) 63 seuset <- Seurat::CreateSeuratObject(counts = counts, min.cells = min_cells, min.features = min_genes)
64 64
65 #+ echo = FALSE 65 #+ echo = FALSE
66 if(showcode == TRUE && vlnfeat == TRUE){print("Raw data vizualization")} 66 if(showcode == TRUE && vlnfeat == TRUE) print("Raw data vizualization")
67 #+ echo = `showcode`, warning = `warn`, include=`vlnfeat` 67 #+ echo = `showcode`, warning = `warn`, include=`vlnfeat`
68 Seurat::VlnPlot(object = seuset, features = c("nFeature_RNA", "nCount_RNA"), axis="v") 68 Seurat::VlnPlot(object = seuset, features = c("nFeature_RNA", "nCount_RNA"))
69 Seurat::FeatureScatter(object = seuset, feature1 = "nCount_RNA", feature2 = "nFeature_RNA") 69 Seurat::FeatureScatter(object = seuset, feature1 = "nCount_RNA", feature2 = "nFeature_RNA")
70 70
71 #+ echo = FALSE 71 #+ echo = FALSE
72 if(showcode == TRUE){print("Filter and normalize for UMI counts")} 72 if(showcode == TRUE) print("Filter and normalize for UMI counts")
73 #+ echo = `showcode`, warning = `warn` 73 #+ echo = `showcode`, warning = `warn`
74 seuset <- subset(seuset, subset = `nCount_RNA` > low_thresholds & `nCount_RNA` < high_thresholds) 74 seuset <- subset(seuset, subset = `nCount_RNA` > low_thresholds & `nCount_RNA` < high_thresholds)
75 seuset <- Seurat::NormalizeData(seuset, normalizeation.method = "LogNormalize", scale.factor = 10000) 75 seuset <- Seurat::NormalizeData(seuset, normalization.method = "LogNormalize", scale.factor = 10000)
76 76
77 #+ echo = FALSE 77 #+ echo = FALSE
78 if(showcode == TRUE && featplot == TRUE){print("Variable Genes")} 78 if(showcode == TRUE && featplot == TRUE) print("Variable Genes")
79 #+ echo = `showcode`, warning = `warn`, include = `featplot` 79 #+ echo = `showcode`, warning = `warn`, include = `featplot`
80 seuset <- Seurat::FindVariableFeatures(object = seuset, selection.method = "mvp") 80 seuset <- Seurat::FindVariableFeatures(object = seuset, selection.method = "mvp")
81 Seurat::VariableFeaturePlot(seuset, cols = c("black", "red"), selection.method = "disp") 81 Seurat::VariableFeaturePlot(seuset, cols = c("black", "red"), selection.method = "disp")
82 seuset <- Seurat::ScaleData(object = seuset, vars.to.regress = "nCount_RNA") 82 seuset <- Seurat::ScaleData(object = seuset, vars.to.regress = "nCount_RNA")
83 83
84 #+ echo = FALSE 84 #+ echo = FALSE
85 if(showcode == TRUE && PCplots == TRUE){print("PCA Visualization")} 85 if(showcode == TRUE && PCplots == TRUE) print("PCA Visualization")
86 #+ echo = `showcode`, warning = `warn`, include = `PCplots` 86 #+ echo = `showcode`, warning = `warn`, include = `PCplots`
87 seuset <- Seurat::RunPCA(seuset, npcs=numPCs) 87 seuset <- Seurat::RunPCA(seuset, npcs = numPCs)
88 Seurat::VizDimLoadings(seuset, dims = 1:2) 88 Seurat::VizDimLoadings(seuset, dims = 1:2)
89 Seurat::DimPlot(seuset, dims = c(1,2), reduction="pca") 89 Seurat::DimPlot(seuset, dims = c(1,2), reduction = "pca")
90 Seurat::DimHeatmap(seuset, dims=1:numPCs, nfeatures=30, reduction="pca") 90 Seurat::DimHeatmap(seuset, dims = 1:numPCs, nfeatures = 30, reduction = "pca")
91 seuset <- Seurat::JackStraw(seuset, dims=numPCs, reduction = "pca", num.replicate = 100) 91 seuset <- Seurat::JackStraw(seuset, dims=numPCs, reduction = "pca", num.replicate = 100)
92 seuset <- Seurat::ScoreJackStraw(seuset, dims = 1:numPCs) 92 seuset <- Seurat::ScoreJackStraw(seuset, dims = 1:numPCs)
93 Seurat::JackStrawPlot(seuset, dims = 1:numPCs) 93 Seurat::JackStrawPlot(seuset, dims = 1:numPCs)
94 Seurat::ElbowPlot(seuset, ndims = numPCs, reduction = "pca") 94 Seurat::ElbowPlot(seuset, ndims = numPCs, reduction = "pca")
95 95
96 #+ echo = FALSE 96 #+ echo = FALSE
97 if(showcode == TRUE && tsne == TRUE){print("tSNE")} 97 if(showcode == TRUE && tsne == TRUE) print("tSNE")
98 #+ echo = `showcode`, warning = `warn`, include = `tsne` 98 #+ echo = `showcode`, warning = `warn`, include = `tsne`
99 seuset <- Seurat::FindNeighbors(object = seuset) 99 seuset <- Seurat::FindNeighbors(object = seuset)
100 seuset <- Seurat::FindClusters(object = seuset) 100 seuset <- Seurat::FindClusters(object = seuset)
101 seuset <- Seurat::RunTSNE(seuset, dims = 1:numPCs, resolution = resolution) 101 seuset <- Seurat::RunTSNE(seuset, dims = 1:numPCs, resolution = resolution)
102 Seurat::DimPlot(seuset, reduction="tsne") 102 Seurat::DimPlot(seuset, reduction = "tsne")
103 103
104 #+ echo = FALSE 104 #+ echo = FALSE
105 if(showcode == TRUE && heatmaps == TRUE){print("Marker Genes")} 105 if(showcode == TRUE && heatmaps == TRUE) print("Marker Genes")
106 #+ echo = `showcode`, warning = `warn`, include = `heatmaps` 106 #+ echo = `showcode`, warning = `warn`, include = `heatmaps`
107 markers <- Seurat::FindAllMarkers(seuset, only.pos = TRUE, min.pct = min_pct, logfc.threshold = logfc_threshold) 107 markers <- Seurat::FindAllMarkers(seuset, only.pos = TRUE, min.pct = min_pct, logfc.threshold = logfc_threshold)
108 top10 <- dplyr::group_by(markers, cluster) 108 top10 <- dplyr::group_by(markers, cluster)
109 top10 <- dplyr::top_n(top10, 10, avg_logFC) 109 top10 <- dplyr::top_n(top10, 10, avg_logFC)
110 Seurat::DoHeatmap(seuset, features = top10$gene) 110 Seurat::DoHeatmap(seuset, features = top10$gene)