BSmooth.tstat issue

[MBeyens]

Hi

I'm testing the bsseq package on an in-house RRBS WGS dataset. During the bsseq::BSmooth.tstat function I encountered some errors. How can I fix this?

R command

bsseq_fit_tstat <- bsseq::BSmooth.tstat(bsseq_fit_filtered,
                                    group1 = c(sampleNames(bsseq_fit)[grep("treatment", sampleNames(bsseq_fit))]), 
                                    group2 = c(sampleNames(bsseq_fit)[grep("control", sampleNames(bsseq_fit))]),
                                    estimate.var = "group2", # "same", "paired'
                                    local.correct = TRUE,
                                    mc.cores = 20,
                                    verbose = TRUE)

R error

[BSmooth.tstat] preprocessing ... done in 3.5 sec
[BSmooth.tstat] computing stats within groups ... done in 5.1 sec
[BSmooth.tstat] computing stats across groups ... done in 2.5 sec
Warning messages:
1: In mclapply(clusterIdx, compute.correction, qSd = qSd, mc.cores = mc.cores) :
  all scheduled cores encountered errors in user code
2: In cbind(...) :
  number of rows of result is not a multiple of vector length (arg 6)

Running on a 94G/20core node.

R sessionInfo

Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 7 (Core)

Matrix products: default
BLAS/LAPACK: /opt/ebsofts/Compiler/GCC/8.2.0-2.31.1/OpenBLAS/0.3.5/lib/libopenblas_sandybridgep-r0.3.5.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

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

other attached packages:
 [1] HDF5Array_1.12.2            rhdf5_2.28.0               
 [3] bsseq_1.20.0                SummarizedExperiment_1.14.1
 [5] DelayedArray_0.10.0         BiocParallel_1.18.1        
 [7] matrixStats_0.54.0          Biobase_2.44.0             
 [9] GenomicRanges_1.36.0        GenomeInfoDb_1.20.0        
[11] IRanges_2.18.2              S4Vectors_0.22.0           
[13] BiocGenerics_0.30.0        

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.1               compiler_3.6.0           XVector_0.24.0          
 [4] R.methodsS3_1.7.1        bitops_1.0-6             R.utils_2.8.0           
 [7] tools_3.6.0              DelayedMatrixStats_1.6.0 zlibbioc_1.30.0         
[10] lattice_0.20-38          BSgenome_1.52.0          Matrix_1.2-17           
[13] GenomeInfoDbData_1.2.1   rtracklayer_1.44.3       Biostrings_2.52.0       
[16] gtools_3.8.1             locfit_1.5-9.1           grid_3.6.0              
[19] data.table_1.12.2        XML_3.98-1.19            limma_3.40.6            
[22] Rhdf5lib_1.6.0           scales_1.0.0             Rsamtools_2.0.0         
[25] GenomicAlignments_1.20.1 permute_0.9-5            colorspace_1.4-1        
[28] RCurl_1.95-4.12          munsell_0.5.0            R.oo_1.22.0 

Many thanks!

Matthias

[PeteHaitch]

It's a little hard to know from that output.
Does the same code work if you restrict mc.cores to a smaller value?
You might test using a subset of bsseq_fit_filtered, e.g, bsseq_fit_filtered[1:1000, ].