Exporting to mzML

[ococrook]

I'm trying to export a spectra object with backend MsBackendMzR() to a mzML file as given by the example in the vignette. However, I get an error which is a bit cryptic. Sorry Johannes I'm no trying to break everying on purpose ....

exampleIonMob1

MSn data (Spectra) with 25000 spectra in a MsBackendMzR backend:
        msLevel     rtime scanIndex
      <integer> <numeric> <integer>
1             1   574.303    307662
2             2   295.926    130990
3             2   659.956    362070
4             1   644.558    352381
5             2   234.457     91995
...         ...       ...       ...
24996         2   214.332     79055
24997         1   612.014    331754
24998         2   199.876     69909
24999         2   429.040    215512
25000         2   159.767     44531
 ... 33 more variables/columns.

file(s):
191115_PhosB_HDMSe_0sec_01_vial01.mzML

fl <- tempfile()
export(exampleIonMob1, MsBackendMzR(), file = fl)

Error in z[1L, 1L] : subscript out of bounds
21.
stop(e)
20.
value[[3L]](cond)
19.
tryCatchOne(expr, names, parentenv, handlers[[1L]])
18.
tryCatchList(expr, classes, parentenv, handlers)
17.
tryCatch({ FUN(...) }, error = handle_error)
16.
withCallingHandlers({ tryCatch({ FUN(...) }, error = handle_error) ...
15.
FUN(...)
14.
FUN(X[[i]], ...)
13.
lapply(X, FUN_, ...)
12.
bplapply(X, FUN, ..., BPREDO = BPREDO, BPPARAM = param)
11.
bplapply(X, FUN, ..., BPREDO = BPREDO, BPPARAM = param)
10.
bplapply(X = seq_along(ddd[[1L]]), wrap, .FUN = FUN, .ddd = ddd, .MoreArgs = MoreArgs, BPREDO = BPREDO, BPPARAM = BPPARAM)
9.
bplapply(X = seq_along(ddd[[1L]]), wrap, .FUN = FUN, .ddd = ddd, .MoreArgs = MoreArgs, BPREDO = BPREDO, BPPARAM = BPPARAM)
8.
bpmapply(.write_ms_data_mzR, split(x, f), levels(f), MoreArgs = list(format = format, copy = copy), BPPARAM = BPPARAM)
7.
bpmapply(.write_ms_data_mzR, split(x, f), levels(f), MoreArgs = list(format = format, copy = copy), BPPARAM = BPPARAM)
6.
.local(object, ...)
5.
export(backend, object, ...)
4.
export(backend, object, ...)
3.
.local(object, ...)
2.
export(exampleIonMob1, MsBackendMzR(), file = fl)
1.
export(exampleIonMob1, MsBackendMzR(), file = fl)

sessionInfo()

 version 4.0.2 (2020-06-22)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 18362)

Matrix products: default

locale:
[1] LC_COLLATE=English_United Kingdom.1252  LC_CTYPE=English_United Kingdom.1252    LC_MONETARY=English_United Kingdom.1252
[4] LC_NUMERIC=C                            LC_TIME=English_United Kingdom.1252    

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

other attached packages:
 [1] forcats_0.5.0        stringr_1.4.0        dplyr_1.0.0          purrr_0.3.4          readr_1.4.0         
 [6] tidyr_1.1.2          tibble_3.0.1         ggplot2_3.3.2        tidyverse_1.3.0      Spectra_0.99.8      
[11] MASS_7.3-51.6        pheatmap_1.0.12      pRoloc_1.29.0        BiocParallel_1.23.2  MLInterfaces_1.69.0 
[16] cluster_2.1.0        annotate_1.67.1      XML_3.99-0.5         AnnotationDbi_1.51.3 IRanges_2.23.10     
[21] MSnbase_2.15.7       ProtGenerics_1.21.0  S4Vectors_0.27.14    mzR_2.23.1           Rcpp_1.0.4.6        
[26] Biobase_2.49.1       BiocGenerics_0.35.4 

loaded via a namespace (and not attached):
  [1] backports_1.1.10      readxl_1.3.1          BiocFileCache_1.13.1  plyr_1.8.6            splines_4.0.2        
  [6] inline_0.3.16         digest_0.6.25         htmltools_0.5.0       foreach_1.5.0         viridis_0.5.1        
 [11] fansi_0.4.1           magrittr_1.5          memoise_1.1.0         doParallel_1.0.15     mixtools_1.2.0       
 [16] openxlsx_4.1.5        limma_3.45.14         recipes_0.1.13        modelr_0.1.8          gower_0.2.2          
 [21] matrixStats_0.57.0    RcppParallel_5.0.2    askpass_1.1           lpSolve_5.6.15        prettyunits_1.1.1    
 [26] colorspace_1.4-1      rvest_0.3.6           blob_1.2.1            rappdirs_0.3.1        haven_2.2.0          
 [31] xfun_0.18             callr_3.5.1           jsonlite_1.7.1        crayon_1.3.4          RCurl_1.98-1.2       
 [36] hexbin_1.28.1         lme4_1.1-23           impute_1.63.0         survival_3.1-12       iterators_1.0.12     
 [41] glue_1.4.1            gtable_0.3.0          ipred_0.9-9           zlibbioc_1.35.0       V8_3.2.0             
 [46] pkgbuild_1.1.0        car_3.0-10            kernlab_0.9-29        rstan_2.21.2          abind_1.4-5          
 [51] scales_1.1.1          vsn_3.57.0            mvtnorm_1.1-1         DBI_1.1.0             viridisLite_0.3.0    
 [56] xtable_1.8-4          progress_1.2.2        foreign_0.8-80        bit_4.0.4             proxy_0.4-24         
 [61] mclust_5.4.6          preprocessCore_1.51.0 StanHeaders_2.21.0-6  MsCoreUtils_1.1.7     lava_1.6.8           
 [66] prodlim_2019.11.13    sampling_2.8          httr_1.4.2            FNN_1.1.3             RColorBrewer_1.1-2   
 [71] ellipsis_0.3.1        loo_2.3.1             farver_2.0.3          pkgconfig_2.0.3       nnet_7.3-14          
 [76] dbplyr_1.4.4          caret_6.0-86          labeling_0.3          tidyselect_1.1.0      rlang_0.4.6          
 [81] reshape2_1.4.4        munsell_0.5.0         cellranger_1.1.0      tools_4.0.2           LaplacesDemon_16.1.4 
 [86] cli_2.1.0             generics_0.0.2        RSQLite_2.2.1         broom_0.7.1           evaluate_0.14        
 [91] mzID_1.27.0           yaml_2.2.1            processx_3.4.4        fs_1.5.0              ModelMetrics_1.2.2.2 
 [96] knitr_1.30            bit64_4.0.5           zip_2.0.4             randomForest_4.6-14   dendextend_1.14.0    
[101] ncdf4_1.17            nlme_3.1-148          xml2_1.3.2            biomaRt_2.45.5        compiler_4.0.2       
[106] rstudioapi_0.11       curl_4.3              e1071_1.7-3           affyio_1.59.0         reprex_0.3.0         
[111] statmod_1.4.34        stringi_1.4.6         ps_1.4.0              lattice_0.20-41       Matrix_1.2-18        
[116] nloptr_1.2.2.1        vctrs_0.3.0           pillar_1.4.6          lifecycle_0.2.0       BiocManager_1.30.10  
[121] MALDIquant_1.19.3     data.table_1.12.8     bitops_1.0-6          R6_2.4.1              pcaMethods_1.81.0    
[126] affy_1.67.1           gridExtra_2.3         rio_0.5.16            codetools_0.2-16      boot_1.3-25          
[131] gtools_3.8.2          assertthat_0.2.1      openssl_1.4.3         withr_2.3.0           hms_0.5.3            
[136] grid_4.0.2            rpart_4.1-15          timeDate_3043.102     coda_0.19-4           class_7.3-17         
[141] minqa_1.2.4           rmarkdown_2.4         segmented_1.2-0       carData_3.0-4         pROC_1.16.2          
[146] lubridate_1.7.9  

[jorainer]

Indeed - pretty cryptic. Could you please disable parallel processing first by setting register(SerialParam()) and then try again. That way we can hopefully see the real error message (I guess something fails in mzR - we never tested writing ion mobility data).

[ococrook]

Unfortunatly no less crytic...

register(SerialParam())
registered()

$SerialParam
class: SerialParam
  bpisup: FALSE; bpnworkers: 1; bptasks: 0; bpjobname: BPJOB
  bplog: FALSE; bpthreshold: INFO; bpstopOnError: TRUE
  bpRNGseed: ; bptimeout: 2592000; bpprogressbar: FALSE
  bpexportglobals: TRUE
  bplogdir: NA
  bpresultdir: NA

$SnowParam
class: SnowParam
  bpisup: FALSE; bpnworkers: 6; bptasks: 0; bpjobname: BPJOB
  bplog: FALSE; bpthreshold: INFO; bpstopOnError: TRUE
  bpRNGseed: ; bptimeout: 2592000; bpprogressbar: FALSE
  bpexportglobals: TRUE
  bplogdir: NA
  bpresultdir: NA
  cluster type: SOCK

fl <- tempfile()
export(exampleIonMob1, MsBackendMzR(), file = fl, BPPARAM = SerialParam())

19.
stop(e)
18.
value[[3L]](cond)
17.
tryCatchOne(expr, names, parentenv, handlers[[1L]])
16.
tryCatchList(expr, classes, parentenv, handlers)
15.
tryCatch({ FUN(...) }, error = handle_error)
14.
withCallingHandlers({ tryCatch({ FUN(...) }, error = handle_error) ...
13.
FUN(...)
12.
FUN(X[[i]], ...)
11.
lapply(X, FUN_, ...)
10.
bplapply(X = seq_along(ddd[[1L]]), wrap, .FUN = FUN, .ddd = ddd, .MoreArgs = MoreArgs, BPREDO = BPREDO, BPPARAM = BPPARAM)
9.
bplapply(X = seq_along(ddd[[1L]]), wrap, .FUN = FUN, .ddd = ddd, .MoreArgs = MoreArgs, BPREDO = BPREDO, BPPARAM = BPPARAM)
8.
bpmapply(.write_ms_data_mzR, split(x, f), levels(f), MoreArgs = list(format = format, copy = copy), BPPARAM = BPPARAM)
7.
bpmapply(.write_ms_data_mzR, split(x, f), levels(f), MoreArgs = list(format = format, copy = copy), BPPARAM = BPPARAM)
6.
.local(object, ...)
5.
export(backend, object, ...)
4.
export(backend, object, ...)
3.
.local(object, ...)
2.
export(exampleIonMob1, MsBackendMzR(), file = fl, BPPARAM = SerialParam())
1.
export(exampleIonMob1, MsBackendMzR(), file = fl, BPPARAM = SerialParam())

[jorainer]

OK, indeed not very helpful - could you eventually provide an example file so that I can reproduce locally?

[ococrook]

I sent an email to johannes.rainer@eurac.edu with a google drive link with 2 ion mobility datasets. Let me know if you got it?

[jorainer]

Did not get any email yet...

[ococrook]

Nothing with the subject line "ion mobility package"?

[jorainer]

Sorry, yes, sure - I found the mail now. And export fails with any of the two files you shared, right?

[ococrook]

No worries, yes neither of them work.

[jorainer]

I finally had time to look into the files. I checked the file "2017121001.mzML". The problem with that is that I've got spectra in there without any peaks. Could you please check this also for your files from above?

sum(lengths(exampleIonMob1) == 0)

If that's the case I have to think how we're going to handle that.

[jorainer]

Note: I've fixed the empty spectra bug in both the BioC 3.12 release branch and the current development (master/main) branch. Could you please install any of the two to see if it fixes your problem @ococrook (e.g. with devtools::install_github("RforMassSpectrometry/Spectra", ref = "RELEASE_3_12"))?

[ococrook]

Great, thank you so much. I'll have a check later

[ococrook]

Hi @jorainer, managed to get around to check and looks good to me! My local versions exported correctly. Indeed, there are empty spectra - about 6025/25000 were empty.

Just to confirm, if we have Hdf5 backend, we first convert to mzR backend to export.

Happy to close the issue.

[jorainer]

If you have a Hdf5 backend you can call export(sps, backend = MsBackendMzR()), so there is no need to change the backend first. With backend you can specify which backend should be used for the export - which will also define the output format.