Mercurial > repos > iuc > scater_plot_pca

comparison scater-create-qcmetric-ready-sce.R @ 0:bea3359ba852 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:11:45 -0400
parents
children 9e5c0bb18d08
comparison
equal deleted inserted replaced
-1:000000000000 0:bea3359ba852
1 #!/usr/bin/env Rscript
2 #Creates a SingleCellExperiment object, which scater's calculateQCMetrics already applied
3
4 library(optparse)
5 library(workflowscriptscommon)
6 library(scater)
7 library(LoomExperiment)
8
9 # parse options
10 #SCE-specific options
11 option_list = list(
12 make_option(
13 c("-a", "--counts"),
14 action = "store",
15 default = NA,
16 type = 'character',
17 help = "A tab-delimited expression matrix. The first column of all files is assumed to be feature names and the first row is assumed to be sample names."
18 ),
19 make_option(
20 c("-r", "--row-data"),
21 action = "store",
22 default = NULL,
23 type = 'character',
24 help = "Path to TSV (tab-delimited) format file describing the features. Row names from the expression matrix (-a), if present, become the row names of the SingleCellExperiment."
25 ),
26 make_option(
27 c("-c", "--col-data"),
28 action = "store",
29 default = NULL,
30 type = 'character',
31 help = "Path to TSV format file describing the samples (annotation). The number of rows (samples) must equal the number of columns in the expression matrix."
32 ),
33 #The scater-specific options
34 make_option(
35 c("--assay-name"),
36 action = "store",
37 default = 'counts',
38 type = 'character',
39 help= "String specifying the name of the 'assay' of the 'object' that should be used to define expression."
40 ),
41 make_option(
42 c("-f", "--mt-controls"),
43 action = "store",
44 default = NULL,
45 type = 'character',
46 help = "Path to file containing a list of the mitochondrial control genes"
47 ),
48 make_option(
49 c("-p", "--ercc-controls"),
50 action = "store",
51 default = NULL,
52 type = 'character',
53 help = "Path to file containing a list of the ERCC controls"
54 ),
55 make_option(
56 c("-l", "--cell-controls"),
57 action = "store",
58 default = NULL,
59 type = 'character',
60 help = "Path to file (one cell per line) to be used to derive a vector of cell (sample) names used to identify cell controls (for example, blank wells or bulk controls)."
61 ),
62 make_option(
63 c("-o", "--output-loom"),
64 action = "store",
65 default = NA,
66 type = 'character',
67 help = "File name in which to store the SingleCellExperiment object in Loom format."
68 )
69 )
70
71 opt <- wsc_parse_args(option_list, mandatory = c('counts', 'output_loom'))
72
73 # Read the expression matrix
74
75 counts <- wsc_split_string(opt$counts)
76 reads <- read.table(counts)
77
78 # Read row and column annotations
79
80 rowdata <- opt$row_data
81
82 if ( ! is.null(opt$row_data) ){
83 rowdata <- read.delim(opt$row_data)
84 }
85
86 coldata <- opt$col_data
87
88 if ( ! is.null(opt$col_data) ){
89 coldata <- read.delim(opt$col_data)
90 }
91
92 # Now build the object
93 assays <- list(as.matrix(reads))
94 names(assays) <- c(opt$assay_name)
95 scle <- SingleCellLoomExperiment(assays = assays, colData = coldata, rowData = rowdata)
96 # Define spikes (if supplied)
97
98
99 #Scater options
100
101 # Check feature_controls (only mitochondrial and ERCC used for now)
102 feature_controls_list = list()
103 if (! is.null(opt$mt_controls) && opt$mt_controls != 'NULL'){
104 if (! file.exists(opt$mt_controls)){
105 stop((paste('Supplied feature_controls file', opt$mt_controls, 'does not exist')))
106 } else {
107 mt_controls <- readLines(opt$mt_controls)
108 feature_controls_list[["MT"]] <- mt_controls
109 }
110 }
111
112 if (! is.null(opt$ercc_controls) && opt$ercc_controls != 'NULL'){
113 if (! file.exists(opt$ercc_controls)){
114 stop((paste('Supplied feature_controls file', opt$ercc_controls, 'does not exist')))
115 } else {
116 ercc_controls <- readLines(opt$ercc_controls)
117 feature_controls_list[["ERCC"]] <- ercc_controls
118 }
119 } else {
120 ercc_controls <- character()
121 }
122
123 # Check cell_controls
124 cell_controls_list <- list()
125 if (! is.null(opt$cell_controls) && opt$cell_controls != 'NULL'){
126 if (! file.exists(opt$cell_controls)){
127 stop((paste('Supplied feature_controls file', opt$cell_controls, 'does not exist')))
128 } else {
129 cell_controls <- readLines(opt$cell_controls)
130 cell_controls_list[["empty"]] <- cell_controls
131 }
132 }
133
134
135 # calculate QCMs
136 scle <- calculateQCMetrics(scle, exprs_values = opt$assay_name, feature_controls = feature_controls_list, cell_controls = cell_controls_list)
137
138 # Output to a Loom file
139 if (file.exists(opt$output_loom)) {
140 file.remove(opt$output_loom)
141 }
142 export(scle, opt$output_loom, format='loom')