CRAM support
Opened this issue · 6 comments
By @GWarsow:
Conversion to cram files requires aligned and sorted sam/bam files.
During the generation of the final sorted/merged/duplicate-marked bam file, it makes no sense to format to cram.
Hence, we stay with bam format for the alignment, sorting, duplicate-tagging and merging (nothing changes in bwaMemSortSlim.sh)
After merging and duplicate-marking, we could think of going to cram.
But as we already have the sam/bam stream available after starting from fastq files, we can use them for downstream analyses.
Important to state: We should avoid writing bam files to HD and write only the cram file instead.
Concludingly, the alignment/QC workflow does not have to switch to cram but only has to create a cram file as last step.
By @GWarsow:
Things may get more complicated when the merged cram file already exists and new lanes shall be added to it.
samtools view will correctly give back singlebams that are already inside the cram file.
In order to contue using the script alreadyMergedLane.pl, the parameter pairedBamSuffix has to be adapted ( e.g. "_paired.bam.sorted.bam" --> "_paired.bam.sorted.cram") or pairedCramSuffix has to be introduced.
By @GWarsow:
Now lets look at marking duplicates. As a first step, we will skip the Picard-based block and concentrate on the biobambam block.
If the cram file already exists, we have to fill the QC NPs with BAM content (and SAM content for NP_COMBINEDANALYSIS_IN.
This can be achieved with cramtools/samtools. But we have to see, whether this works out well with regard to performance. Obviously, just performing a cat on the bam file will not suffice when we have a cram file.
First implementation is like
cat ${FILENAME} | ${CRAM2BAM_TARGET} | tee ${NP_FLAGSTATS_IN}...
See https://eilslabs-phabricator.dkfz.de/T280 for more details.
CRAM generation should use samtools, which is lossless by default. The way to go for completely lossless production and restoring CRAM is then probably
# cram->bam
samtools view -O cram,store_md=1,store_nm=1 -o aln.cram aln.bam
# bam->cram
samtools view --input-fmt cram,decode_md=0 -o aln.new.bam aln.cram This will even store and restore the NM and MD tags.
An index will have to be build separately, like with bam by
samtools index aln.cram aln.crai (see samtools man-page).