Review and optimise the number of file conversions

[muffato]

Description of feature

The pipeline is undergoing a major overhaul with the introduction of chunking.
We need to review the number of conversions it does (between BAM, CRAM, and FASTQ) and try to minimise those. This is to save space and runtime.

[muffato]

#211 and #184 did an excellent job for the alignment part.

The remaining unnecessary conversions are:

• CRAM to FASTQ for FastQC because it doesn't support CRAM
• In CONVERT_STATS: running crumble, converting to CRAM, changing the header. All create a new file.

[muffato]

Process list for a test_full run:

INPUT_CHECK:SAMTOOLS_FLAGSTAT (ilPolIcar4.ont.fastq)                                                       4 of 4 ✔
PREPARE_GENOME:GUNZIP (GCA_937595015.1.fasta.gz)                                                           1 of 1 ✔
PREPARE_GENOME:MASK_UNMASK (GCA_937595015.1)                                                               1 of 1 ✔
PREPARE_GENOME:BWAMEM2_INDEX (GCA_937595015.1.unmasked.fa)                                                 1 of 1 ✔
FASTQC (SAMEA5928167.RO_PI_748)                                                                            4 of 4 ✔

ALIGN_SHORT:CONVERT_CRAM                                                                                   -
ALIGN_SHORT:SAMTOOLS_ADDREPLACERG                                                                          -
ALIGN_SHORT:CRAM_MAP_ILLUMINA:SAMTOOLS_INDEX                                                               -
ALIGN_SHORT:CRAM_MAP_ILLUMINA:CRAMALIGN_GENCRAMCHUNKS (2)                                                  2 of 2 ✔
ALIGN_SHORT:CRAM_MAP_ILLUMINA:SAMTOOLS_SPLITHEADER (SAMEA7523244.ERR9123833)                               2 of 2 ✔
ALIGN_SHORT:CRAM_MAP_ILLUMINA:BWAMEM2_INDEX (GCA_937595015.1.unmasked.fa)                                  2 of 2 ✔
ALIGN_SHORT:CRAM_MAP_ILLUMINA:CRAMALIGN_BWAMEM2ALIGNHIC (SAMEA7523244.ERR9123832)                          5 of 5 ✔
ALIGN_SHORT:CRAM_MAP_ILLUMINA:BAM_SAMTOOLS_MERGE_MARKDUP:SAMTOOLS_FAIDX (GCA_937595015.1.unmasked.fa)      2 of 2 ✔
ALIGN_SHORT:CRAM_MAP_ILLUMINA:BAM_SAMTOOLS_MERGE_MARKDUP:SAMTOOLS_MERGEDUP (SAMEA7523244.ERR9123833)       2 of 2 ✔

ALIGN_LONG:SAMTOOLS_SPLITHEADER (SAMEA7523243.ERR9081703)                                                  1 of 1 ✔
ALIGN_LONG:GAWK_MODIFY_YAML_BARCODE (SAMEA7523243.ERR9081703)                                              1 of 1 ✔
ALIGN_LONG:PACBIO_PREPROCESS_ULI:LIMA                                                                      -
ALIGN_LONG:PACBIO_PREPROCESS_ULI:PBMARKDUP                                                                 -
ALIGN_LONG:PACBIO_PREPROCESS_ULI:BLAST_MAKEBLASTDB (HiFi_adapter)                                          1 of 1 ✔
ALIGN_LONG:PACBIO_PREPROCESS_ULI:BLAST_BLASTN                                                              -
ALIGN_LONG:PACBIO_PREPROCESS_ULI:HIFITRIMMER_PROCESSBLAST                                                  -
ALIGN_LONG:PACBIO_PREPROCESS_ULI:HIFITRIMMER_FILTERBAM                                                     -
ALIGN_LONG:PACBIO_PREPROCESS:BLAST_MAKEBLASTDB (HiFi_adapter)                                              1 of 1 ✔
ALIGN_LONG:PACBIO_PREPROCESS:BLAST_BLASTN (SAMEA7523243.ERR9081703)                                        1 of 1 ✔
ALIGN_LONG:PACBIO_PREPROCESS:HIFITRIMMER_PROCESSBLAST (SAMEA7523243.ERR9081703)                            1 of 1 ✔
ALIGN_LONG:PACBIO_PREPROCESS:HIFITRIMMER_FILTERBAM (SAMEA7523243.ERR9081703)                               1 of 1 ✔
ALIGN_LONG:FASTQC_FILTERED (SAMEA7523243.ERR9081703)                                                       1 of 1 ✔
ALIGN_LONG:SAMTOOLS_FASTQ                                                                                  -
ALIGN_LONG:FASTX_MAP_LONG_READS:FASTXALIGN_PYFASTXINDEX (SAMEA7523243.ERR9081703)                          2 of 2 ✔
ALIGN_LONG:FASTX_MAP_LONG_READS:MINIMAP2_INDEX (2)                                                         2 of 2 ✔
ALIGN_LONG:FASTX_MAP_LONG_READS:FASTXALIGN_MINIMAP2ALIGN (SAMEA7523243.ERR9081703)                         823 of 823 ✔
ALIGN_LONG:FASTX_MAP_LONG_READS:BAM_SAMTOOLS_MERGE_MARKDUP:SAMTOOLS_FAIDX (GCA_937595015.1.unmasked.fa)    2 of 2 ✔
ALIGN_LONG:FASTX_MAP_LONG_READS:BAM_SAMTOOLS_MERGE_MARKDUP:SAMTOOLS_MERGE (SAMEA5928167.RO_PI_748)         2 of 2 ✔

CONVERT_STATS:CRUMBLE (SAMEA7523243.ERR9081703)                                                            3 of 3 ✔
CONVERT_STATS:CHANGE_NAME (SAMEA7523243.ERR9081703)                                                        4 of 4 ✔
CONVERT_STATS:SAMTOOLS_CRAM (SAMEA7523243.ERR9081703)                                                      4 of 4, retries: 2
CONVERT_STATS:BLOBTK_DEPTH (SAMEA7523243.ERR9081703)                                                       4 of 4, retries: 1
CONVERT_STATS:BGZIP_BEDGRAPH (GCA_937595015.1.unmasked.pacbio.SAMEA7523243.ERR9081703.minimap2.regions.bed.gz) 4 of 4 ✔
CONVERT_STATS:SAMTOOLS_STATS (SAMEA5928167.RO_PI_748)                                                      4 of 4 ✔
CONVERT_STATS:GZIP_STATS (SAMEA5928167.RO_PI_748)                                                          4 of 4 ✔
CONVERT_STATS:SAMTOOLS_FLAGSTAT (SAMEA7523243.ERR9081703)                                                  4 of 4 ✔
CONVERT_STATS:SAMTOOLS_IDXSTATS (SAMEA7523243.ERR9081703)                                                  4 of 4 ✔
MULTIQC                                                                                                    1 of 1 ✔

and sizes (in MB)

                IN    OUT  x
<genome>       495  11252 25
SAMEA7523243 17084  72301  5
SAMEA7523244 23249 163181  7
SAMEA5928167 39650 166430  5

Current path for short reads:

• Conversion of the input to CRAM if needed [usually not needed for us]
• CRAM aligned with bwa-mem2, by chunks, stored as BAM (CRAMALIGN_BWAMEM2ALIGNHIC)
• Aligned chunks are merged, still using BAM (SAMTOOLS_MERGEDUP)

Current path for PacBio:

• if hifi trimmer: Removes and filters reads, into FastQ (HIFITRIMMER_FILTERBAM). On the fly of the input from BAM/CRAM if needed [usually the case for us]
• else: convert to FastQ if needed [usually the case for us]
• Aligns the FastQ with minimap2 by chunks, stored as BAM (FASTXALIGN_MINIMAP2ALIGN)
• Aligned chunks are merged, still using BAM (SAMTOOLS_MERGE)

Final processing:

• Compress scores with Crumble, still using BAM (CRUMBLE)
• Convert to CRAM (SAMTOOLS_CRAM)

This results in the following disk space usage (test_full profile):

• 25x the genome size (in Gbp ~ the uncompressed Fasta file).
• 5x the size of the input PacBio files
• 7x the size of the input Illumina/HiC files
• 5x the size of the input ONT files

My thoughts:

1. This could be optimised by using CRAM as the output format of samtools merge-ing the chunks. That would require flipping the bam/cram outputs in CONVERT_STATS. We also need to find where to set the CRAM output options (v3.0, archive mode, etc).
2. Then, we could use CRAM instead of BAM as the intermediate format within the sanger-tol alignment sub-workflows
3. Finally, we could use CRAM as the input format of the long read minimap alignment sub-workflow (and output of the hifi-trimmer sub-workflow)

For intermediate CRAM usage, we'd have to pass the reference Fasta (doable since it's available across the sub-workflows). I recommend using the default CRAM settings as it's already a >2x size reduction over BAM.

cc @sainsachiko