--- a/scater-pca-filter.R	Tue Sep 03 14:30:21 2019 -0400
+++ b/scater-pca-filter.R	Thu Sep 09 12:24:39 2021 +0000
@@ -12,49 +12,48 @@
 library(workflowscriptscommon)
 library(LoomExperiment)
 library(scater)
-library(mvoutlier)
+library(robustbase)
 
 # parse options
-option_list = list(
+option_list <- list(
   make_option(
     c("-i", "--input-loom"),
     action = "store",
     default = NA,
-    type = 'character',
+    type = "character",
     help = "A SingleCellExperiment object file in Loom format."
   ),
   make_option(
     c("-o", "--output-loom"),
     action = "store",
     default = NA,
-    type = 'character',
+    type = "character",
     help = "File name in which to store the SingleCellExperiment object in Loom format."
   )
 )
 
-opt <- wsc_parse_args(option_list, mandatory = c('input_loom', 'output_loom'))
+opt <- wsc_parse_args(option_list, mandatory = c("input_loom", "output_loom"))
 
 # Check parameter values
 
-if ( ! file.exists(opt$input_loom)){
-  stop((paste('File', opt$input_loom, 'does not exist')))
+if (! file.exists(opt$input_loom)) {
+  stop((paste("File", opt$input_loom, "does not exist")))
 }
 
-# Input from Loom format
+# Filter out unexpressed features
 
-scle <- import(opt$input_loom, format='loom', type='SingleCellLoomExperiment')
-print(paste("Starting with", ncol(scle), "cells and", nrow(scle), "features."))
+sce <- import(opt$input_loom, format = "loom", type = "SingleCellLoomExperiment")
+
+print(paste("Starting with", ncol(sce), "cells"))
 
-# Run PCA on data and detect outliers
-scle <- runPCA(scle, use_coldata = TRUE, detect_outliers = TRUE)
+sce <- runColDataPCA(sce, outliers = TRUE, variables = list("sum", "detected", "subsets_Mito_percent"))
+sce$use <- !sce$outlier
+sce <- sce[, colData(sce)$use]
+print(paste("Ending with", ncol(sce), "cells"))
 
-# Filter out outliers
-scle <- scle[, !scle$outlier]
-
-print(paste("Ending with", ncol(scle), "cells and", nrow(scle), "features."))
 
 # Output to a Loom file
 if (file.exists(opt$output_loom)) {
   file.remove(opt$output_loom)
 }
-export(scle, opt$output_loom, format='loom')
+export(sce, opt$output_loom, format = "loom")