Mercurial > repos > iuc > scater_filter

annotate scater-plot-dist-scatter.R @ 0:e6ca62ac65c6 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
author iuc
date Thu, 18 Jul 2019 11:13:41 -0400
parents
children b7ea9f09c02f
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
1 #!/usr/bin/env Rscript
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
2
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
3 # Plot the distribution of read counts and feature counts, side by side, then a scatter plot of read counts vs feature counts below
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
4
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
5 # Load optparse we need to check inputs
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
6
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
7 library(optparse)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
8 library(workflowscriptscommon)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
9 library(LoomExperiment)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
10 library(scater)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
11 library(ggpubr)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
12
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
13 # parse options
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
14
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
15 option_list = list(
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
16 make_option(
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
17 c("-i", "--input-loom"),
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
18 action = "store",
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
19 default = NA,
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
20 type = 'character',
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
21 help = "A SingleCellExperiment object file in Loom format."
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
22 ),
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
23 make_option(
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
24 c("-o", "--output-plot-file"),
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
25 action = "store",
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
26 default = NA,
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
27 type = 'character',
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
28 help = "Path of the PDF output file to save plot to."
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
29 )
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
30 )
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
31
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
32 opt <- wsc_parse_args(option_list, mandatory = c('input_loom', 'output_plot_file'))
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
33
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
34 # Check parameter values
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
35
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
36 if ( ! file.exists(opt$input_loom)){
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
37 stop((paste('File', opt$input_loom, 'does not exist')))
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
38 }
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
39
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
40 # Input from Loom format
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
41
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
42 scle <- import(opt$input_loom, format='loom', type='SingleCellLoomExperiment')
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
43
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
44 #do the scatter plot of reads vs genes
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
45 total_counts <- scle$total_counts
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
46 total_features <- scle$total_features_by_counts
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
47 count_feats <- cbind(total_counts, total_features)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
48 cf_dm <- as.data.frame(count_feats)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
49
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
50 # Calculate binwidths for reads and features plots. Use 20 bins
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
51 read_bins <- max(total_counts / 1e6) / 20
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
52 feat_bins <- max(total_features) / 20
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
53
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
54 #make the plots
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
55 plot <- ggplot(cf_dm, aes(x=total_counts / 1e6, y=total_features)) + geom_point(shape=1) + geom_smooth() + xlab("Read count (millions)") +
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
56 ylab("Feature count") + ggtitle("Scatterplot of reads vs features")
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
57 plot1 <- qplot(total_counts / 1e6, geom="histogram", binwidth = read_bins, ylab="Number of cells", xlab = "Read counts (millions)", fill=I("darkseagreen3")) + ggtitle("Read counts per cell")
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
58 plot2 <- qplot(total_features, geom="histogram", binwidth = feat_bins, ylab="Number of cells", xlab = "Feature counts", fill=I("darkseagreen3")) + ggtitle("Feature counts per cell")
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
59 plot3 <- plotColData(scle, y="pct_counts_MT", x="total_features_by_counts") + ggtitle("% MT genes") + geom_point(shape=1) + theme(text = element_text(size=15)) + theme(plot.title = element_text(size=15))
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
60
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
61 final_plot <- ggarrange(plot1, plot2, plot, plot3, ncol=2, nrow=2)
e6ca62ac65c6 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/scater commit 5fdcafccb6c645d301db040dfeed693d7b6b4278
iuc
parents:
diff changeset
62 ggsave(opt$output_plot_file, final_plot, device="pdf")