Multicore Bismark corrupts result files if one of child process fails

[iromeo]

Hi, I'm using Bismark (0.23.0) in Snakemake pipeline to process many of the WGBS. I've noticed that some of my BAMs have fewer reads than expected, and ambiguous/unaligned FASTQ files have invalid format. I run my bismark in multicore mode. When I checked logs, I noticed that 2 of 4 cores failed with Child process terminated with exit signal: '13' but the main Bismark process continued working and returned 0 exit code. As a result I got an incorrect result which hard to notice in batch processing mode. In this case Bismark should fail with an error. Child cores failed likely due to cluster issues (e.g. lack of HDD space or etc), I successfully re-aligned one of the problematic samples and got correct results and files.

Some details:

• Command to launch Bismark: bismark -X 600 --gzip --multicore 4 --ambiguous --unmapped --bowtie2 $(dirname resources/indexes/hg38/Bisulfite_Genome) -1 data/reads/wgbs_y23_pooled/Clean/FE04/FE04_1.fq.gz -2 data/reads/wgbs_y23_pooled/Clean/FE04/FE04_2.fq.gz
• Interesting pars of log:

Current dir: /scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/shadow/tmp60xczd7h
Bowtie 2 seems to be working fine (tested command 'bowtie2 --version' [2.4.2])
Output format is BAM (default)
Alignments will be written out in BAM format. Samtools found here: '/scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/conda/af71104ea71ba452675d70f85d913ba3/bin/samtools'
Reference genome folder provided is resources/indexes/hg38/     (absolute path is '/storage1/fs1/martyomov/Active/projects/aging/LongitudinalYear123/indexes/hg38/)'
FastQ format assumed (by default)

Input files to be analysed (in current folder '/scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/shadow/tmp60xczd7h'):
data/reads/wgbs_y23_pooled/Clean/FE04/FE04_1.fq.gz
data/reads/wgbs_y23_pooled/Clean/FE04/FE04_2.fq.gz
Library is assumed to be strand-specific (directional), alignments to strands complementary to the original top or bottom strands will be ignored (i.e. not performed!)
Summary of all aligner options: -q --score-min L,0,-0.2 --ignore-quals --no-mixed --no-discordant --dovetail --maxins 600
Running Bismark Parallel version. Number of parallel instances to be spawned: 4

Current working directory is: /scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/shadow/tmp60xczd7h

(...)

C methylated in CHH context:    0.5%
C methylated in unknown context (CN or CHN):    1.0%



Child process terminated with exit signal: '13'


Child process terminated with exit signal: '13'



Right, cleaning up now...

Deleting temporary sequence files...
FE04_1.fq.gz.temp.1.gz  FE04_2.fq.gz.temp.1.gz  FE04_1.fq.gz.temp.2.gz  FE04_2.fq.gz.temp.2.gz  FE04_1.fq.gz.temp.3.gz  FE04_2.fq.gz.temp.3.gz  FE04_1.fq.gz.temp.4.gz  FE04_2.fq.gz.temp.4.gz

Now merging BAM files FE04_1.fq.gz.temp.1.gz_bismark_bt2_pe.bam FE04_1.fq.gz.temp.2.gz_bismark_bt2_pe.bam FE04_1.fq.gz.temp.3.gz_bismark_bt2_pe.bam FE04_1.fq.gz.temp.4.gz_bismark_bt2_pe.bam into >>> FE04_1_bismark_bt2_pe.bam <<<
Merging from file >> FE04_1.fq.gz.temp.1.gz_bismark_bt2_pe.bam <<
Merging from file >> FE04_1.fq.gz.temp.2.gz_bismark_bt2_pe.bam <<
[W::bam_hdr_read] EOF marker is absent. The input is probably truncated
[E::bgzf_read] Read block operation failed with error 4 after 0 of 4 bytes
[main_samview] truncated file.
Failed to close filehandle
Merging from file >> FE04_1.fq.gz.temp.3.gz_bismark_bt2_pe.bam <<
[W::bam_hdr_read] EOF marker is absent. The input is probably truncated
[E::bgzf_read] Read block operation failed with error 4 after 0 of 4 bytes

(...)

Merging from file >> FE04_1.fq.gz.temp.1.gz_bismark_bt2_PE_report.txt <<
Merging from file >> FE04_1.fq.gz.temp.2.gz_bismark_bt2_PE_report.txt <<
Use of uninitialized value $total_seqs in addition (+) at /scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/conda/af71104ea71ba452675d70f85d913ba3/bin/bismark line 1128, <ALN> line 45.
Use of uninitialized value $unique in addition (+) at /scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/conda/af71104ea71ba452675d70f85d913ba3/bin/bismark line 1129, <ALN> line 45.
Use of uninitialized value $no_aln in addition (+) at /scratch1/fs1/martyomov/rcherniatchik/aging-wgbs-smk-pipeline-v3/.snakemake/conda/af71104ea71ba452675d70f85d913ba3/bin/bismark line 1130, <ALN> line 45.

(...)

Deleting temporary report files...
FE04_1.fq.gz.temp.1.gz_bismark_bt2_PE_report.txt        FE04_1.fq.gz.temp.2.gz_bismark_bt2_PE_report.txt        FE04_1.fq.gz.temp.3.gz_bismark_bt2_PE_report.txt        FE04_1.fq.gz.temp.4.gz_bismark_bt2_PE_report.txt

Bismark completed in 3d 5h 45m 29s

====================
Bismark run complete
====================

As a result:

1. I get "601143749" sequence pairs analyzed instead of "1202287499" (in FE04_hg38_PE_report.txt)
2. FASTQ file format is invalid:

Number of lines in *.fq.gz files cannot be divided by 4:

# actual:
results/hg38/bams_bismark/FE04_/FE04_hg38_ambiguous_reads_1.fq.gz - 233298342
results/hg38/bams_bismark/FE04_/FE04_hg38_ambiguous_reads_2.fq.gz - 233298342

# should be:
results/hg38/bams_bismark/FE04/FE04_hg38_ambiguous_reads_1.fq.gz - 277355436
results/hg38/bams_bismark/FE04/FE04_hg38_ambiguous_reads_2.fq.gz - 277355436

Here is the place where format is corrupted:

---

I noticed this file format issue only due to FASTQC error on ambiguous reads. Otherwise the truncated BAM file could go to downstream:

Started analysis of FE04_hg38_ambiguous_reads_1.fq.gz
Approx 5% complete for FE04_hg38_ambiguous_reads_1.fq.gz
Approx 10% complete for FE04_hg38_ambiguous_reads_1.fq.gz
(....)
Approx 40% complete for FE04_hg38_ambiguous_reads_1.fq.gz
Approx 45% complete for FE04_hg38_ambiguous_reads_1.fq.gz
Failed to process file FE04_hg38_ambiguous_reads_1.fq.gz
uk.ac.babraham.FastQC.Sequence.SequenceFormatException: ID line didn't start with '@'
        at uk.ac.babraham.FastQC.Sequence.FastQFile.readNext(FastQFile.java:158)
        at uk.ac.babraham.FastQC.Sequence.FastQFile.next(FastQFile.java:125)
        at uk.ac.babraham.FastQC.Analysis.AnalysisRunner.run(AnalysisRunner.java:77)
        at java.base/java.lang.Thread.run(Thread.java:834)

[iromeo]

Log file
FE04_hg38_pe.bam.log.zip

[FelixKrueger]

Hi Roman,

thanks for reporting this unwanted behaviour in such meticulous detail. It is true that the parallel processing does little (if any) to ensure that the output hasn't been corrupted, e.g. as a result of insufficient disk resources etc.

In that sense, I think it is a good idea to simply kill the job to alert the user that something has gone wrong and needs to be looked at. I shall accept your PR now, and close this issue.