cyttools/logicleTransform.R at master · bc2zb/cyttools · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
#!/usr/bin/env Rscript

require(docopt)
require(methods)

"
Usage:
inspectr.R (-h | --help | --version)
inspectr.R DIR

Description:   This script perfroms QC analysis on spectral flow data, calulating similarity scores and spillover spreading
Options:
--version       Show the current version.

Arguments:

DIR    Provide directory where cyttools.args.Rdata file is located
" -> doc


args <- docopt(doc)

ARGS_DIR <- args$DIR

cat("\nLoading arguments from", ARGS_DIR, "\n")

load(paste(ARGS_DIR, "cyttools.args.Rdata", sep = ""))

RESULTS_DIR <- args$OUT

source("cyttoolsFunctions.R")

##########################################################################
############################ R code goes here ############################
##########################################################################


# capture all files in the directory
all_fcs_files <- list.files(args$DIR, full.names = T, pattern = "\\.fcs$")

# read in fcs files
ncfs <- read.ncdfFlowSet(all_fcs_files)

# transform flowSet
chnls <- colnames(ncfs)[grep("SSC|FSC|Time|\\-H", colnames(ncfs), invert = T)]
safe_estimate_logicle <- safely(estimateLogicle)
transFuncts <- fsApply(ncfs, safe_estimate_logicle, channels = paste0("^", chnls)) %>%
  modify_depth(1, 1) %>%
  discard(is_null)

safe_transform <- safely(transform)

for ( i in 1:length(transFuncts)){
  ncfs_trans <- safe_transform(ncfs, transFuncts[[i]])
  if(is.null(ncfs_trans$error)){
    ncfs_trans <- ncfs_trans$result
    break
  }else if(i == length(transFuncts)){
    cat("\nERROR: No transform can be estimated, exiting now\n")
    q()
  }
}

# preprocess flowFrames

if(args$clean == T){
  ncfs_trans <- fsApply(ncfs_trans, safe_preprocess_frame) %>%
    modify_depth(1, 1) %>%
    discard(is_null) %>%
    flowSet()
}

dir.create(paste0(RESULTS_DIR, "TRANSFORMED_FCS/"),
           showWarnings = F)

for( files in all_fcs_files){
  out.fcs.file <- paste0(RESULTS_DIR, "TRANSFORMED_FCS/logicle_transformed_", basename(files))
  out_fcs_frame <- ncfs_trans[basename(files)][[1]]
  write.FCS(out_fcs_frame, out.fcs.file)
}
##########################################################################
############################     End code     ############################
##########################################################################