Error executing dnmtools pmd (core dumped) probably related to GSL

[ryansohny]

Describe the bug
Hi, thank you for creating this amazing tool.

I'm trying to call Partially Methylated Domain (PMD) using dnmtools with DNA methylation called from PacBio HiFi Reads from which I made input file suitable for calling PMD. Below is the exact command I excuted:

$dnmtools pmd \
-itr 1000 \
-bin 1000 \
-seed 42 \
-verbose \
-o ./PMD/${sample}.pmd \
${sample}_haplotagged_wIndels_splt_mincov4.combined.denovo.ML.meth

and the input file:

head ${sample}_haplotagged_wIndels_splt_mincov4.combined.denovo.ML.meth 
#chr1	10469	+	CpG	75.0 	4
#chr1	10471	+	CpG	75.0 	4
#chr1	10484	+	CpG	100.0 	4
#chr1	10489	+	CpG	100.0 	5
#chr1	10493	+	CpG	100.0 	5
#chr1	10497	+	CpG	100.0 	5
#chr1	10525	+	CpG	80.0 	5
#chr1	10542	+	CpG	100.0 	5
#chr1	10563	+	CpG	100.0 	5
#chr1	10571	+	CpG	100.0 	5

However, I got this error message with a core dumped file:

[DYNAMICALLY SELECTING BIN SIZE]
[READING IN AT BIN SIZE 1000]
[READING CPGS AND METH PROPS] from ${sample}_haplotagged_wIndels_splt_mincov4.combined.denovo.ML.meth
TOTAL BINS: 3019070
MEAN COVERAGE: 281.768
[SEPARATING BY CPG DESERT]
BINS RETAINED: 2876173
NUMBER OF DISTANCES BETWEEN: 2876173
DESERTS REMOVED: 1130

MAX_ITER=1000 TOLERANCE=1e-05
ITR F size B size DELTA
gsl: gamma.c:1561: ERROR: domain error
Default GSL error handler invoked.
/mnt/mone/ge2011.11/bdcm/spool/bdcm05/job_scripts/504341: line 23: 50603 Aborted (core dumped) $dnmtools pmd -itr 1000 -bin 1000 -seed 42 -verbose -o ./PMD/${sample}.pmd ${sample}_haplotagged_wIndels_splt_mincov4.combined.denovo.ML.meth

I first thought the error was brought up because the file had too many lines in it (29066896 number of CpGs), so instead I picked only CpGs in the chromosome 1 (2380932 CpGs) and executed the same command seeing the same error as above.

OS version and total memory:

cat /etc/os-release 
#NAME="CentOS Linux"
#VERSION="7 (Core)"
#ID="centos"
#ID_LIKE="rhel fedora"
#VERSION_ID="7"
#PRETTY_NAME="CentOS Linux 7 (Core)"
#ANSI_COLOR="0;31"
#CPE_NAME="cpe:/o:centos:centos:7"
#HOME_URL="https://www.centos.org/"
#BUG_REPORT_URL="https://bugs.centos.org/"

#CENTOS_MANTISBT_PROJECT="CentOS-7"
#CENTOS_MANTISBT_PROJECT_VERSION="7"
#REDHAT_SUPPORT_PRODUCT="centos"
#REDHAT_SUPPORT_PRODUCT_VERSION="7"

cat /proc/meminfo | grep 'MemTotal'
#MemTotal:       247359496 kB

I also used methpipe (v5.0.1) pmd, and had the same error.

Can you help me with this matter?

Kind Regards,
Sohn

[andrewdavidsmith]

@ryansohny The counts format we use has methylation as fractional. I think you somehow made it into percentages.

https://dnmtools.readthedocs.io/en/latest/counts/

If this turns out to be the only source of your problem, I'll relabel this issue as a need for input validation.

[ryansohny]

Wow this is embarassing. A classic example of "why you should read the manual carefully".
When I created an input from the same sample sequenced in a different way (WGBS), I formatted it in the right way (fractional) so I guess I must have been distracted for some reason.
And I got what I want! It went smoothly.

Thank you very much Professor Smith. Wish you all the best.
Sohn

[andrewdavidsmith]

@ryansohny I appreciate that you brought this up. Today we added a check for this in each place where that file format is loaded (6 different commands I think). The error messages aren't great yet, but at least users will be informed of a problem if those values are not in the expected range.

[ryansohny]

That's great news! Always grateful for developers who go the extra mile to make things easier for their users. Much appreciated it!