methylKit icon indicating copy to clipboard operation
methylKit copied to clipboard

Issues with processBismarkAln()

Open miyakokodama opened this issue 4 years ago • 7 comments

miyakokodama commented 3 minutes ago Hi,

I am new to methylKit, I am currently trying to read my bam files with processBismarkAln however I get an error message I don't understand.

setwd("/home/projects/5_bismark")

file.list=list("BM_174/BM_174.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"BM_427/BM_427.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"BM_437/BM_437.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"BM_455/BM_455.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"E_220/E_220.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"E_256/E_256.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"E_265/E_265.pair1.truncated.gz_bismark_bt2_pe_sorted.bam",
"E_274/E_274.pair1.truncated.gz_bismark_bt2_pe_sorted.bam")

CpG.bismark.data <- processBismarkAln(location=file.list, sample.id = list("BM1","BM2","BM3","BM4","E1","E2","E3","E4"), treatment =c(0,0,0,0,1,1,1,1), assembly = "GCF_000233375.1_ICSASG_v2_genomic", save.folder=getwd(), save.context = c("CpG"), read.context="CpG", nolap=FALSE, mincov=10, minqual=20)

Once these are typed in, I get a following message.


/home/projects/5_bismark/BM_174/BM_174.pair1.truncated.gz_bismark_bt2_pe_sorted.bam
 using htslib.

Error in methCall(read1 = location, type = "bam", nolap = nolap, minqual = minqual,  :
  basic_string::_M_construct null not valid

Any hints/advice you could provide would be greatly appreciated.

Miyako

miyakokodama avatar May 01 '20 15:05 miyakokodama

Hi Miyako,

Could you please send your session info? Also please verify that all your bamfiles are fine, maybe with samtools quickcheck .

Best, Alex

alexg9010 avatar May 02 '20 22:05 alexg9010

Please also provide a small reproducible example, because otherwise I cannot test for your error.

alexg9010 avatar May 03 '20 12:05 alexg9010

Hi Alex, thanks for the reply! I will get back to you as soon as I can.

miyakokodama avatar May 11 '20 10:05 miyakokodama

Hi, Here is the sessioninfo():

> sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 7 (Core)

Matrix products: default
BLAS/LAPACK: /services/tools/intel/perflibs/2019_update5/compilers_and_libraries_2019.5.281/linux/mkl/lib/intel64_lin/libmkl_gf_lp64.so

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

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

other attached packages:
 [1] methylKit_1.12.0     GenomicRanges_1.38.0 GenomeInfoDb_1.22.0
 [4] IRanges_2.20.2       S4Vectors_0.24.3     BiocGenerics_0.32.0
 [7] RCurl_1.95-4.12      bitops_1.0-6         devtools_2.2.2
[10] usethis_1.5.1        sessioninfo_1.1.1    memoise_1.1.0
[13] BiocManager_1.30.10  RLinuxModules_0.2

loaded via a namespace (and not attached):
 [1] Biobase_2.46.0              pkgload_1.0.2
 [3] splines_3.6.1               R.utils_2.9.2
 [5] gtools_3.8.1                assertthat_0.2.1
 [7] GenomeInfoDbData_1.2.2      Rsamtools_2.2.3
 [9] remotes_2.1.1               numDeriv_2016.8-1.1
[11] pillar_1.4.4                backports_1.1.6
[13] lattice_0.20-40             glue_1.4.0
[15] limma_3.42.2                bbmle_1.0.23.1
[17] digest_0.6.25               XVector_0.26.0
[19] qvalue_2.18.0               colorspace_1.4-1
[21] Matrix_1.2-18               R.oo_1.23.0
[23] plyr_1.8.6                  XML_3.98-1.20
[25] pkgconfig_2.0.3             emdbook_1.3.12
[27] zlibbioc_1.32.0             purrr_0.3.4
[29] mvtnorm_1.1-0               scales_1.1.0
[31] processx_3.4.2              BiocParallel_1.20.1
[33] tibble_3.0.1                ggplot2_3.3.0
[35] ellipsis_0.3.0              SummarizedExperiment_1.16.1
[37] withr_2.2.0                 cli_2.0.2
[39] mclust_5.4.5                magrittr_1.5
[41] crayon_1.3.4                ps_1.3.2
[43] R.methodsS3_1.8.0           fs_1.3.2
[45] fansi_0.4.1                 MASS_7.3-51.5
[47] pkgbuild_1.0.7              tools_3.6.1
[49] data.table_1.12.8           prettyunits_1.1.1
[51] matrixStats_0.56.0          lifecycle_0.2.0
[53] stringr_1.4.0               munsell_0.5.0
[55] DelayedArray_0.12.2         callr_3.4.3
[57] Biostrings_2.54.0           compiler_3.6.1
[59] fastseg_1.32.0              rlang_0.4.6
[61] grid_3.6.1                  testthat_2.3.2
[63] gtable_0.3.0                reshape2_1.4.3
[65] R6_2.4.1                    GenomicAlignments_1.22.1
[67] rtracklayer_1.46.0          dplyr_0.8.5
[69] bdsmatrix_1.3-4             rprojroot_1.3-2
[71] desc_1.2.0                  stringi_1.4.6
[73] Rcpp_1.0.4.6                vctrs_0.2.4
[75] tidyselect_1.0.0            coda_0.19-3

I also did samtools quickcheck but no errors were detected. However, when I used picard ValidateSamFile, I got the following error:

WARNING 2020-05-13 13:20:33     ValidateSamFile NM validation cannot be performed without the reference. All other validations will still occur.
ERROR: Read groups is empty
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:1090:1000_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:1090:1000_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:24306:1031_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:24306:1031_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:17526:1094_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:17526:1094_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:9182:1110_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:9182:1110_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:26738:1110_1:N:0:TGTCTAGA+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:26738:1110_1:N:0:TGTCTAGA+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:22064:1125_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:22064:1125_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:3251:1141_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:3251:1141_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:3369:1157_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:3369:1157_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:28492:1172_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:28492:1172_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:11993:1188_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:11993:1188_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:6090:1204_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:6090:1204_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:26756:1204_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:26756:1204_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:4065:1235_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:4065:1235_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:8901:1282_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:8901:1282_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:27525:1282_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:27525:1282_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:27615:1313_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:27615:1313_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:21603:1329_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:21603:1329_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:15519:1344_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:15519:1344_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:7392:1360_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:7392:1360_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:24985:1360_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:24985:1360_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:24722:1376_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:24722:1376_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:14353:1391_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:14353:1391_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:2031:1407_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:2031:1407_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:22263:1438_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:22263:1438_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:13774:1454_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:13774:1454_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:10583:1470_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:10583:1470_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:12373:1501_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:12373:1501_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:19036:1517_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:19036:1517_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:2736:1532_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:2736:1532_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:11677:1548_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:11677:1548_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:31367:1548_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:31367:1548_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:23610:1579_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:23610:1579_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:1542:1626_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:1542:1626_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:13892:1626_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:13892:1626_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:31286:1626_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:31286:1626_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:6451:1705_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:6451:1705_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:14742:1720_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:14742:1720_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:25283:1752_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:25283:1752_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:21296:1767_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:21296:1767_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:3766:1783_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:3766:1783_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:13883:1799_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:13883:1799_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:22354:1814_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:22354:1814_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:12599:1830_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:12599:1830_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:14543:1846_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:14543:1846_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:15980:1861_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:15980:1861_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:27028:1861_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:27028:1861_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:1316:1892_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:1316:1892_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:15447:1908_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:15447:1908_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:16206:1939_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:16206:1939_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
WARNING: Read name A01045:96:HWNN7DSXX:4:1101:31286:1939_1:N:0:TGTCTAGC+GCTACTCA, A record is missing a read group
Maximum output of [100] errors reached.

I would also like to add that processBismarkAln() produces an output ("CpG.txt" and "CpG_conversionStats.txt" files) for the first sample, and then quits with the error I mentioned in the first post. The "CpG.txt" seems to be truncated.

If any of these errors don't ring a bell, please let me know - I will try to subset my .bam files and provide a small reproducible example.

Many thanks for your help in advance!

Miyako

miyakokodama avatar May 13 '20 12:05 miyakokodama

Hi @miyakokodama ,

please provide me a subset of bam files that causes the error, because this does not seem to be an obvious error. I cannot say wether this is a bug in the code or wether something is wrong with your samples unless I have the small reproducible example.

Best, Alex

alexg9010 avatar May 19 '20 08:05 alexg9010

Thanks Alex! I will do that and get back to you.

miyakokodama avatar May 19 '20 19:05 miyakokodama

Hi, dear @alexg9010 @miyakokodama , Now I also suffer this issue with processBismarkAln(), so could you please share the solution with me again?

guobin8216 avatar Sep 27 '22 21:09 guobin8216