Segmentation fault when processing reads

[amleckenby]

Hey, I am getting a segmentation fault when using the latest version of graph map ( e8d6100 )

I am invoking graph map using the command:

./bin/Linux-x64/graphmap align -C -r EhR2_v2.fasta -d all_subreads.fastq -o graphmap_out.sam

I have ran the command many times and find that the program aborts on differing reads whilst processing reads (i.e. sometimes it will segfault at read 0/141336 and what appears at random numbers between read 0 to around read 40/141336).

Here is my output.

[15:26:29 Index] Running in normal (parsimonious) mode. Only one index will be used. [15:26:29 Index] Index already exists. Loading from file. [15:26:30 Index] Index loaded in 0.17 sec. [15:26:30 Index] Memory consumption: [currentRSS = 259 MB, peakRSS = 259 MB]

[15:26:30 Run] Automatically setting the maximum allowed number of regions: max. 500, attempt to reduce after 0 [15:26:30 Run] No limit to the maximum number of seed hits will be set in region selection. [15:26:30 Run] Reference genome is assumed to be circular. [15:26:30 Run] Only one alignment will be reported per mapped read. [15:26:30 ProcessReads] Reads will be loaded in batches of up to 1024 MB in size. [15:26:35 ProcessReads] Batch of 141336 reads (1024 MiB) loaded in 5.22 sec. (12151080 bases) [15:26:35 ProcessReads] Memory consumption: [currentRSS = 1334 MB, peakRSS = 1334 MB] [15:26:35 ProcessReads] Using 24 threads. [15:26:35 ProcessReads] [CPU time: 5.24 sec, RSS: 1334 MB] Read: 8/141336 (0.01%) [m: 0, u: 0], length = 8250, qname: m150916_104438_42215_c10085833255000000182...Segmentation fault

The last line of output differs each time (i.e. a different read is quoted) the program is invoked and subsequently segfaults.

Could you possibly help with a solution?

NB: our reads are pacbio reads.

[isovic]

Hi, thank you for the detailed report!

Would it be possible that you share a few of your reads (and the reference) so I can try to reproduce this issue?

The reason why the read id is different every time is multithreading - the terminal displays only what's on the first thread, and which read gets assigned to which thread is a race condition.

You can try re-running the process in one thread with -t 1, in this case when program segfaults the qname shown in terminal should be the last one that was supposed to have been processed but failed. Could you extract this read into a separate file, and try mapping only this read to check whether it breaks as expected? If it does, could you send it to me for further inspection? Alternatively, could you extract the first 100-500 reads from your dataset and share those? The segfault should occur in those based on your information.

Thank you, Best regards, Ivan.

[amleckenby]

Hi Ivan,

Thanks for your suggestions. I mapped it using -t 1 and it is the same read causing the segfault (75th in the file). I then mapped this read individually (segfaultread.fastq) with no problems so I am still confused as to why the read continuously segfaults in the larger fast file.

I have attached the reference and first 150 reads for you. I have also included the read file which causes the segfault.

Thanks, Amber

[isovic]

Hi, I did some long-overdue profiling in the meantime, and detected a memory access bug. This is probably related to your segfault. I'll merge the update soon (as soon as some other updates are ready as well) and let you know when it's ready.

Best regards, Ivan.

[isovic]

Hi Amber, could you please try pulling the newest commit and testing if it still segfaults?

Thank you, Best regards, Ivan.

[amleckenby]

Hi Ivan,

Tried rerunning the command and now the program terminates after throwing an instance of 'std::bad::alloc'.

See output: ./bin/Linux-x64/graphmap align -C -r ~amber/CANU_assembly/rDNA_analysis/EhR2_v2.fasta -d ~amber/CANU_assembly/all_subreads.fasta -o all_subreads2EhR2_v2.sam

[14:37:20 Index] Running in normal (parsimonious) mode. Only one index will be used. [14:37:20 Index] Index already exists. Loading from file. [14:37:20 Index] Index loaded in 0.29 sec. [14:37:20 Index] Memory consumption: [currentRSS = 259 MB, peakRSS = 259 MB]

[14:37:20 Run] Automatically setting the maximum allowed number of regions: max. 500, attempt to reduce after 0 [14:37:20 Run] No limit to the maximum number of seed hits will be set in region selection. [14:37:20 Run] Reference genome is assumed to be circular. [14:37:20 Run] Only one alignment will be reported per mapped read. [14:37:20 ProcessReads] Reads will be loaded in batches of up to 1024 MB in size. [14:37:27 ProcessReads] Batch of 282026 reads (1024 MiB) loaded in 7.29 sec. (37641976 bases) [14:37:27 ProcessReads] Memory consumption: [currentRSS = 1380 MB, peakRSS = 1380 MB] [14:37:27 ProcessReads] Using 24 threads. [14:37:27 ProcessReads] [CPU time: 7.31 sec, RSS: 1380 MB] Read: 13/282026 (0.00%)[14:37:29 ProcessReads] [CPU time: 59.56 sec, RSS: 1382 MB] Read: 29/282026 (0.01%) [m: 2, u: 4], length = 3819, qname: m150916_104438_42215_c100858332550000001...terminate called after throwing an instance of 'std::bad_alloc' what(): std::bad_alloc Aborted

I am confused by this, I am using a very large computer, normally used for multiple denovo genome assemblies in parallel so the 1TB of RAM should be more than enough for this program?

Amber

[isovic]

Hi Amber, I'm sorry that it still segfaults. I'd like inspect it myself, but I noticed your earlier message was missing the attachment you mentioned. Could you by any chance resend the attachments? P.S. GitHub might have a problem with certain extensions in attachmets - if it bugs you, just rename the extension to .txt instead of .fastq. Thank you, Best regards, Ivan.

[amleckenby]

Hey Ivan,

Sure. Ive attached the reference and the reads as .txt. The reference should be a .fasta and the reads are .fastq

Amber

first_150_subreads.txt ref.txt

[isovic]

Got them now! :-) Will inspect and let you know.

Best regards, Ivan.

[isovic]

Hi Amber, Indeed there were a few memory issues in v0.3.2 which were especially evident on your data. These should now be resolved, I profiled the newest version and checked with your data, and it seems to work fine now. Could you try pulling the newest version (v0.4.1) and verifying that it works?

Thank you! Best regards, Ivan.

[amleckenby]

Hi Ivan,

I can confirm no change in the program performance.

Still aborting due to throwing std::bad_alloc

Program seems to run fine if the -C flag is dropped from the command so my guess would be something is throwing the program by telling it the reference is circular.

output: align -C -r EhR2_v2.fasta -d all_subreads.fastq -o all_subreads_2_EhR2v2_graphmap.sam -t 1 [13:08:11 Index] Running in normal (parsimonious) mode. Only one index will be used. [13:08:11 BuildIndex] Index already exists. Loading from file. [13:08:11 Index] Index loaded in 0.23 sec. [13:08:11 Index] Memory consumption: [currentRSS = 259 MB, peakRSS = 259 MB]

[13:08:11 Run] Automatically setting the maximum allowed number of regions: max. 500, attempt to reduce after 0 [13:08:11 Run] No limit to the maximum number of seed hits will be set in region selection. [13:08:11 Run] Reference genome is assumed to be circular. [13:08:11 Run] Only one alignment will be reported per mapped read. [13:08:11 ProcessReads] Reads will be loaded in batches of up to 1024 MB in size. [13:08:16 ProcessReads] Batch of 141336 reads (1024 MiB) loaded in 4.89 sec. (38784824 bases) [13:08:16 ProcessReads] Memory consumption: [currentRSS = 1335 MB, peakRSS = 1335 MB] [13:08:16 ProcessReads] Using 1 threads. [13:08:16 ProcessReads] [CPU time: 4.89 sec, RSS: 1335 MB] Read: 0/141336 (0.00%) [m: 0, u: 0], length = 5054, qname: m150916_104438_42215_c10085833255000000182...terminate called after throwing an instance of 'std::bad_alloc' what(): std::bad_alloc Aborted amber@gauss04:~/graphmap$