# HG changeset patch # User iuc # Date 1535315042 14400 # Node ID 8d8412d35247a7815c32da5c7340dcbc54b7e0a8 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/seurat commit 24c0223b9baa6d59bba381ef94f7e77b1c204d80 diff -r 000000000000 -r 8d8412d35247 seurat.R --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/seurat.R Sun Aug 26 16:24:02 2018 -0400 @@ -0,0 +1,116 @@ +options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } ) + +# we need that to not crash galaxy with an UTF8 error on German LC settings. +loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") + +suppressPackageStartupMessages({ + library(Seurat) + library(SingleCellExperiment) + library(dplyr) + library(optparse) +}) + +option_list <- list( + make_option(c("-counts","--counts"), type="character", help="Counts file"), + make_option(c("-numPCs","--numPCs"), type="integer", help="Number of PCs to use in plots"), + make_option(c("-min.cells","--min.cells"), type="integer", help="Minimum cells to include"), + make_option(c("-min.genes","--min.genes"), type="integer", help="Minimum genes to include"), + make_option(c("-low.thresholds","--low.thresholds"), type="double", help="Low threshold for filtering cells"), + make_option(c("-high.thresholds","--high.thresholds"), type="double", help="High threshold for filtering cells"), + make_option(c("-x.low.cutoff","--x.low.cutoff"), type="double", help="X-axis low cutoff for variable genes"), + make_option(c("-x.high.cutoff","--x.high.cutoff"), type="double", help="X-axis high cutoff for variable genes"), + make_option(c("-y.cutoff","--y.cutoff"), type="double", help="Y-axis cutoff for variable genes"), + make_option(c("-cells.use","--cells.use"), type="integer", help="Cells to use for PCHeatmap"), + make_option(c("-resolution","--resolution"), type="double", help="Resolution in FindClusters"), + make_option(c("-min.pct","--min.pct"), type="double", help="Minimum percent cells in FindClusters"), + make_option(c("-logfc.threshold","--logfc.threshold"), type="double", help="LogFC threshold in FindClusters"), + make_option(c("-rds","--rds"), type="logical", help="Output Seurat RDS object") + ) + +parser <- OptionParser(usage = "%prog [options] file", option_list=option_list) +args = parse_args(parser) + +counts <- read.delim(args$counts, row.names=1) +seuset <- CreateSeuratObject(raw.data = counts, min.cells = args$min.cells, min.genes = args$min.cells) + +# Open PDF for plots +pdf("out.pdf") + +VlnPlot(object = seuset, features.plot = c("nGene", "nUMI"), nCol = 2) +GenePlot(object = seuset, gene1 = "nUMI", gene2 = "nGene") + +print("Filtering cells") +if (!is.null(args$low.thresholds)){ + lowthresh <- args$low.thresholds +} else { + lowthresh <- "-Inf" +} +if (!is.null(args$high.thresholds)){ + highthresh <- args$high.thresholds +} else { + highthresh <- "Inf" +} +seuset <- FilterCells(object = seuset, subset.names = c("nUMI"), + low.thresholds=c(lowthresh), high.thresholds = c(highthresh)) + +print("Normalizing the data") +seuset <- NormalizeData(object = seuset, normalization.method = "LogNormalize", + scale.factor = 10000) + +print("Finding variable genes") +seuset <- FindVariableGenes(object = seuset, mean.function = ExpMean, + dispersion.function = LogVMR, + x.low.cutoff = args$x.low.cutoff, + x.high.cutoff = args$x.high.cutoff,, + y.cutoff = args$y.cutoff +) + +print("Scaling the data and removing unwanted sources of variation") +seuset <- ScaleData(object = seuset, vars.to.regress = c("nUMI")) + +print("Performing PCA analysis") +seuset <- RunPCA(object = seuset, pc.genes = seuset@var.genes) +VizPCA(object = seuset, pcs.use = 1:2) +PCAPlot(object = seuset, dim.1 = 1, dim.2 = 2) +PCHeatmap( + object = seuset, + pc.use = 1:args$numPCs, + cells.use = args$cell.use, + do.balanced = TRUE, + label.columns = FALSE, + use.full = FALSE +) + +print("Determining statistically significant principal components") +seuset <- JackStraw(object = seuset, num.replicate = 100, display.progress= FALSE) +JackStrawPlot(object = seuset, PCs = 1:args$numPCs) +PCElbowPlot(object = seuset) + +print("Clustering the cells") +seuset <- FindClusters( + object = seuset, + reduction.type = "pca", + dims.use = 1:args$numPCs, + resolution = args$resolution, + print.output = 0, + save.SNN = TRUE +) + +print("Running non-linear dimensional reduction (tSNE)") +seuset <- RunTSNE(object = seuset, dims.use = 1:args$numPCs, do.fast = TRUE) +TSNEPlot(object = seuset) + +print("Finding differentially expressed genes (cluster biomarkers)") +markers <- FindAllMarkers(object = seuset, only.pos = TRUE, min.pct = args$min.pct, + logfc.threshold = args$logfc.threshold) +top10 <- markers %>% group_by(cluster) %>% top_n(10, avg_logFC) +DoHeatmap(object = seuset, genes.use = top10$gene, slim.col.label = TRUE, remove.key = TRUE) + +# Close PDF for plots +dev.off() + +if (!is.null(args$rds) ) { + saveRDS(seuset, "Seurat.rds") +} + +sessionInfo() \ No newline at end of file diff -r 000000000000 -r 8d8412d35247 seurat.xml --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/seurat.xml Sun Aug 26 16:24:02 2018 -0400 @@ -0,0 +1,146 @@ + + - toolkit for exploration of single-cell RNA-seq data + + r-base + r-seurat + bioconductor-singlecellexperiment + r-optparse + + /dev/null | grep -v -i "WARNING: ")", SingleCellExperiment version" $(R --vanilla --slave -e "library(SingleCellExperiment); cat(sessionInfo()\$otherPkgs\$SingleCellExperiment\$Version)" 2> /dev/null | grep -v -i "WARNING: ")", optparse version" $(R --vanilla --slave -e "library(optparse); cat(sessionInfo()\$otherPkgs\$optparse\$Version)" 2> /dev/null | grep -v -i "WARNING: ") + ]]> + + + + + + + +
+ + + + + + + + + + + + +
+
+ + + + + rscript + + + rds + + + + + + + + + + + + + + + + + + + + + + + + + 10.1038/nbt.4096 + +
diff -r 000000000000 -r 8d8412d35247 test-data/deng_small.tab.gz Binary file test-data/deng_small.tab.gz has changed diff -r 000000000000 -r 8d8412d35247 test-data/out.pdf Binary file test-data/out.pdf has changed