ce_get_transformations mixes up scales between channels

[bjreisman]

Hi RGLab,

@sierrabarone and I came across the following bug in the ce_get_transformations function where it appears transformations are not being defined correctly. I think it has something to do with how the co-factor is defined for each channel and whether or not the transformation is used before the function moves on to the next channel. I was able to get it to work as desired in a bit of hacky way, but there is probably a more elegant solution. I've included an reproducible example below (acs file in a zip at the bottom). I'd be happy to create a PR with my hacky fix (calling each transformation once) if you'd like.

-Ben

> library(flowCore)
> library(CytoML)
> library(dplyr)
> library(flowWorkspace)
As part of improvements to flowWorkspace, some behavior of
GatingSet objects has changed. For details, please read the section
titled "The cytoframe and cytoset classes" in the package vignette:

  vignette("flowWorkspace-Introduction", "flowWorkspace")
> 
> list.files() #find the name of the acs file here
[1] "20210408_scales_repex.R"                   "experiment_36838_Apr-08-2021_07-54-PM.acs"
[3] "get_transformations_repex_acs.zip"        
> ce <- CytoML::open_cytobank_experiment("experiment_36838_Apr-08-2021_07-54-PM.acs")
Unpacking ACS file...
> transformers <- ce_get_transformations(ce)
> 
> scales <- bind_rows(lapply(ce$experiment$scales, as.data.frame))
> scales
    channelShortName cofactor maximum minimum scaleType
1              fsc-a        1  262144       1    Linear
2              fsc-h        1  262144       1    Linear
3              fsc-w        1  262144       1    Linear
4              ssc-a        1  262144       1    Linear
5              ssc-h        1  262144       1    Linear
6              ssc-w        1  262144       1    Linear
7  alexa fluor 647-a      300   45000       1   Arcsinh
8  alexa fluor 700-a     3000  262144   -4000   Arcsinh
9           apc-h7-a      150  262144     -10   Arcsinh
10 alexa fluor 488-a    20000  262144     900   Arcsinh
11     percp-cy5-5-a      150  262144   -3000   Arcsinh
12    pacific blue-a      150  200000      10   Arcsinh
13        qdot 525-a      150  262144    -200   Arcsinh
14  pacific orange-a      150  262144    -200   Arcsinh
15        qdot 605-a      150  262144    -200   Arcsinh
16        qdot 655-a      150  262144    -200   Arcsinh
17        qdot 705-a      150  262144    -200   Arcsinh
18              pe-a      150  262144    -200   Arcsinh
19    pe-texas red-a      150  262144    -200   Arcsinh
20           7-aad-a      150  262144    -200   Arcsinh
21    alexa 680-pe-a      150  262144    -200   Arcsinh
22          pe-cy7-a     2000  262144   -2000   Arcsinh
23              time        1  262144       1    Linear
> 
> #should be ~0.88 not 0.075
> transformers$`Pacific Blue-A`$transform(150)
[1] 0.07492986
> 
> #seems to think the co-factor is 2000?
> transformers$`Pacific Blue-A`$transform(2000)
[1] 0.8813736
> 
> 
> 
> # Alternative get_transformations-----------------------------------------------
> ce_get_transformations2 <- function (x) 
+ {
+   chnls <- colnames(x)
+   low.chnls <- tolower(chnls)
+   scales <- x$experiment$scales
+   res <- list()
+   for (scale in scales) {
+     pname <- scale$channelShortName
+     ind <- match(pname, low.chnls)
+     if (is.na(ind)) 
+       stop(pname, " not found in colnames of the panel!")
+     pname <- chnls[ind]
+     stype <- scale$scaleType
+     cofactor <- scale$cofactor
+     if (stype == "Log") 
+       res[[pname]] <- logtGml2_trans()
+     else if (stype == "Arcsinh") {
+       T <- sinh(1) * cofactor
+       res[[pname]] <- asinhtGml2_trans(T = T, M = 0.434294481903252, 
+                                        A = 0)
+       test <- res[[pname]]$transform(cofactor) # this is the line I added
+     }
+     else if (stype != "Linear") 
+       stop("Unknown scale type: ", stype)
+   }
+   transformerList(names(res), res)
+ }
> 
> transformers2 <- ce_get_transformations2(ce)
> 
> #now it seems to work correctly
> transformers2$`Pacific Blue-A`$transform(150)
[1] 0.8813736

> sessionInfo()
R version 4.0.4 (2021-02-15)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19042)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252   
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] flowWorkspace_4.3.2 dplyr_1.0.4         CytoML_2.2.1        flowCore_2.2.0     

loaded via a namespace (and not attached):
 [1] tidyselect_1.1.0    purrr_0.3.4         lattice_0.20-41     colorspace_2.0-0    vctrs_0.3.6        
 [6] generics_0.1.0      stats4_4.0.4        yaml_2.2.1          ncdfFlow_2.36.0     base64enc_0.1-3    
[11] utf8_1.1.4          RBGL_1.66.0         XML_3.99-0.5        rlang_0.4.10        hexbin_1.28.2      
[16] pillar_1.5.0        glue_1.4.2          DBI_1.1.1           aws.s3_0.3.21       Rgraphviz_2.34.0   
[21] BiocGenerics_0.36.0 RColorBrewer_1.1-2  plyr_1.8.6          matrixStats_0.58.0  jpeg_0.1-8.1       
[26] lifecycle_1.0.0     zlibbioc_1.36.0     RProtoBufLib_2.3.1  munsell_0.5.0       gtable_0.3.0       
[31] cytolib_2.3.2       latticeExtra_0.6-29 Biobase_2.50.0      parallel_4.0.4      curl_4.3           
[36] fansi_0.4.2         Rcpp_1.0.6          scales_1.1.1        S4Vectors_0.28.1    jsonlite_1.7.2     
[41] RcppParallel_5.0.2  graph_1.68.0        gridExtra_2.3       ggplot2_3.3.3       png_0.1-7          
[46] digest_0.6.27       grid_4.0.4          tools_4.0.4         magrittr_2.0.1      tibble_3.0.6       
[51] crayon_1.4.1        aws.signature_0.6.0 pkgconfig_2.0.3     ellipsis_0.3.1      data.table_1.14.0  
[56] xml2_1.3.2          assertthat_0.2.1    httr_1.4.2          R6_2.5.0            ggcyto_1.18.0      
[61] compiler_4.0.4   

get_transformations_repex_acs.zip

[mikejiang]

Thanks for reporting this!I've pushed a more permanent fix to flowWorkspace.
see RGLab/flowWorkspace@c29fce7
now all should be good

> t <- sinh(1) * 150
> trans <- asinhtGml2_trans(T = t, M = 0.434294481903252,A = 0)
> t <- sinh(1) * 2000
> trans$transform(150)
[1] 0.8813736